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ABSTRACT 

This study was concerned with the dynamic flexural response 
of propeller blades subjected to harmonic forces, and presents 
both analytical and experimental results. 

The determination of the torsional constants of three blade 
models having NACA four -digit symmetrical airfoil cross sections 
is presented. Values were obtained for these models analytically 
and experimentally. In addition, results were ohtained for three 
other models having rectangular, elliptical, and parabolic cross 
sections. 

Complete modal analyses were performed for five blade models, 
The identification of modal parameters was done for cases when the 
blades were modeled as either undamped or damped multi-degree-of- 
freedom systems. For the experimental phase of this study, the 
modal testing was performed using a Dual Channel FFT analyzer 
and an impact hammer (which produced an impulsive excitation).. The 
natural frequency and damping of each mode in the frequency range 
up to two kilohertz were measured, 

A small computer code was developed to calculate the dynamic 
response of the blade models for comparison with the experimental 
results. A comparison of the undamped and damped cases was 
made for all five blade models at the instant of maximum excitation 
force. The program was capable of handling models where the 
excitation forces were distributed arbitrarily along the length 


of the blade. 
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Chapte/.* 1 
INTRODUCTION 

1.1. General Remarks and Objectives of the Study 

The purpose of this investigation was to study, both analytically 
and experimentally, the dynamics of cantilevered beams having 
airfoil-like cross sections, including damping and inertial 
coupling. Information from this study will be used in a more 
comprehensive analysis of twisted, rotating propeller blades. 

The vibration and fatigue of propeller blades are important from 
both a comfort and safety viewpoint. 

In a study of the dynamics of propeller blades, one must 
consider several factors, including: nonuniform mass and torsional 

inertia distributions; nonuniform flapwxse, chordwise and torsioual 
stiffness distributions; built-in twist; noncoincidence of the elastic 
axis with the body axis; and inertial coupling. Thus, before 
considering the complete case of propeller blade dynamics, it is 
instructive to study the behavior of simpler blade models. 

Specifically, the objectives of the present study were; 

1. To develop a method to determine the torsional stiffness 

of blades having airfoil cross sections. First, an analytical model 
was developed based on existing energy methods in solid mechanics. 
Secondly, these methods were applied to three models having NACA 
four-digit symmetrical airfoil sections, and were verified by 
experimental results. 

2. To define analytically the dynamic flexural response of 
blade models, modeled as either an undamped or damped multi-degree- 


of-freedom system. The amplitudes of harmonic excitation force 
were arbitrarily distributed along the length of the blade. 
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3. To undertake an analytical investigation of dynamic 
characteristics of mult i-degree-of -freedom systems in a frequency 
domain, and their relation to the time domain. Special attention 
was given to the identification of damping. 

4. To determine experimentally the natural frequencies and 
modal damping of models discussed above with a modern FFT analyzer 
using an impulsive excitation. 

5 . To develop a small computer code to calculate either an 
undamped or damped response of blade models discussed above. 

Special attention was paid to the modeling of damping, 

6. To evaluate the results obtained by the computer program. 
Special emphasis was placed on damping and its influence on the 
dynamics of frequencies close to the natural frequency. 

1.2. Previous Investigations on Torsional Stiffness of 
Models Having Airfoil Cross Sections 

Since the torsion of models having a cross section different 
than that of a circle is governed by Poisson's partial differential 
equation, the closed form solution can be obtained only for a few 
regular cases. For configurations not defined analytically, torsional 
problems can be solved only by approximate methods such as the 
finite element and finite difference methods. 

One of the objectives of this study was to develop analytical 
and experimental procedures for determining torsional constants for 
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beams having airfoil-like cross sections. In researching this 
problem, the existing literature was studied. 

1,2.1. Reference Books. Papers, and Literature Studied 

A very elaborate treatment of energy methods in stress analysis 
is given by Richards in reference 62. 

Dym and Shames (14) and Crandal (11) present rigorous treatments 
of torsional problems by variational approaches. Emphasis is placed 
on the development of basic concepts of energy methods for torsional 
problems presented and solutions for models having rectangular, 
quadratic, and parabolic segment cross sections. 

Chou and Pagano (7.) , Kartog (19) , and Oder, and Ripperberg 
(55) present detailed analogies of problems related to torsion. The 
analogies provide physical insight into problems of this nature 
and suggest experimental approaches to these problems. 

Most of the analytical work done In this present study was 
based on the above references. 

Timoshenko and Goodier (86) present a very elaborate method 
for determining the total energy of beams in torsion having an 
arbitrary cross section. They then apply the Rayleigh-Ritz method 
to problems having various cross sections. Basically, the total 
energy is found In terms of stress functions. These stress functions 
are then expanded by a series In order to satisfy the boundary 
conditions. Finally, the extremum for total energy is found with 
respect to parameters in the series. These researchers also present 
this technique for the beams having an airfoil cross section. 


§ 
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Cook (9) gives an extensive treatment of the method of weighted 
residuals for a finite element approach. Since this present study 
does not use a finite element approach, only Cook's concepts were used. 

Finlayson and Scriven (15) and Hutton and Anderson (26) present 
the method of weighted residuals. They give a constructive 
critique on several methods, with special attention given to Galerkln's 
method, 

Duncan (13) shows that the analytical approach to the Saint- 
Venant torsion problem can be readily done by using Galerkin's 
analytical method. He presents a detailed consideration of torsion 
of beams having a symmetrical airfoil and a symmetrical parabolic 
segment as cross sections. Furthermore, he suggests that this 
approach can be employed for nonsymmetrical cross sections by 
using a doubly Infinite series of the function representing the 
contour of the airfoil. Although presented in 1938, this method 
still appears useful for studying the torsion of propeller and rotor 
blades. 

The references discussed above were used as the basis for the 
numerical calculation of torsional stiffness for the three models 
already discussed. 

1.3. Reference Books. Papers, and Literature Studied 
On the Dynamical Response of Propeller Blades 

Kenedy and Pancy (31) present one of the first expositions on 
determining modal characteristics from test data. They assume that 
damping takes a rather specialized form, and use polar plots to 
determine the modal characteristics. 
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Lewis and Wrisley (42) present one of the first multiple-shaker 
systems for modal testing, They use 24 Independently controlled 
shakers to (hopefully) drive the system into a pure normal mode. 

They assume that the shakers are distributed with respect to the 
mass, and that the damping is distributed in the same manner as the 
mass. They also assume that the damping is equal in all modes of 
vibration. 

Trail and Nash (87) propose a multiple-excitation technique 
which makes use of Kenedy and Pancy's studies to locate the natural 
frequencies and the sets of linearly independent forces applied to 
the structure. This is done in order to calculate the required 
forces needed to excite a pure normal mode. 

Asher (1) applies essentially the same analysis as Trail and 
Nash (87) , but suggests an Improved method of locating the natural 
frequencies. He uses the determinant of the real part of the 
flexibility matrix to define the natural frequencies. 

Walgrave and Ehlbeck (90) present a good review of modal 
analysis which gives results in the time and frequency domains. 

Caughey (6), Clough and Penzien (8), and Paz (56) present the 
analytical modeling of modal damping which decouples the damping 
forces. 

Klosterman (32, 33, 34, 35, and 36) gives a very detailed 
treatment on various subjects in analytical and experimental modal 
analysis. Some of his papers are more oriented to experimental 
analysis. 

Craig (10), Meirovitch (49), Bathe (3), and Paz (56) give 
rigorous consideration to various subjects on structural dynamics, 
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Including both classical and numerical approaches. They are oriented 
toward treatment of vibration phenomena in the time domain. Most of 
them were used during the consideration of the undamped and damped 
multi-degree-of -freedom systems in Chapters 5 and 6 of this study. 

Richardson (63, 64, 65, 66, 67, and 68), Potter (58,59) and 
Ramsey (60) cover the treatment of vibrational systems having 
multiple degrees of freedom. The authors present concepts of 
vibration analysis in the frequency domain, and they are pioneering 
efforts in treating the vibration of multi-degree-of -freedom systems 
via .system dynamics. 

In order to Improve accuracy with modern FFT analyzers, the 
following steps should be taken: 

1. Calibration of the impact hammer. 

2. Conditioning of the input and output signals. 

3. Consideration of the noise and its Influence on the 
measured transfer function. 

4. Proper selection of the frequency range so that all 
frequencies are excited within the chosen range. 

The Spectral Dynamics Co. (82), Harris (17,18), the Wavetek-Rockland 
Co. (91), the Nicolet Scientific Co. (54), and the Hewlett-Packard 
Co. (21) present detailed treatments of signals in. order to properly 
conduct vibration testing. 

An important point in an experimental program is the accurate 
measurement of the structural transfer function. Here noise effects 
must be considered. Keller (28,29,30), Lally (37,38), Brown (5), 
Mitchel and Lynch (50) , and Halvorsen and Bendat (16) present 
various subjects of crucial interest for determining structural 



transfer functions. They give a detailed treatment of noise in 
the time and frequency domains and its influence during transfer 
function measurements. Also, they discuss the calibration of an 
impact hammer, and input and output signals. 

Finally, measured data must be Interpreted properly. 

Lang (39,40), the Spectral Dynamics Co. (80,81), the Nicolet 
Scientific Co. (51, 52, 53), the Hewlett -Packard Co, (22), and 
Brown and Halvorsen (5) present in detail the Interpretation of 
data obtained by an FFT analyzer. 


Chapter 2 


ANALYTICAL METHODS FOR DETERMINATION OF TORSIONAL 
STIFFNESS .OF BEAMS HAVING AIRFOIL-LIKE CROSS SECTIONS 

2,1, Torsional Behavior of Beams Having an 
Arbitrary Cross Section 

In the case of the torsion of beams not circular in cross 
section, we lose the arguments of symmetry and, along with them, 
the simplicity of the elementary theory of torsion. The argument 
that plane cross sections remain plane during deformation, for 
example, is now no longer valid. Observation shows that the cross 
section of a non-circular section does not remain plane upon 
twisting, but warps out of its plane. 

Other than the fact that cross sections warp, perhaps one of 

the most obvious characteristics of the beam’s behavior is the 

absence of normal stresses. No external forces or bending moments 

are present, and no end constraints exist; therefore, the only 

stress components needed to provide the equilibrium Of any 

transverse segment are shearing stress ‘in the cross-sectional 

planes. Furthermore, of the three components of shearing stress, 

only T z y and can result in a twisting moment: The component 

T is zero. Thus, we conclude that for the noncircular beam in 
xy 

torsion. 



The differential equations of equilibrium in terms of stresses 
and body forces reduce to 
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( 2 , 2 ) 

(2.3) 


Equation (2.3) shows that the shearing stresses do not vary with 
z and, hence, have the same distribution on each cross section. 
Observing eq. (2.2), we have a statically Indeterminate problem, 
and we are forced to turn to considerations of strains and dis- 
placements for additional Information. Since the material is 
assumed to be homogenous Isotropic, and linearly elastic, we may 
Introduce eq.. (2.1) Into equations for stress and strain in a 
homogenous isotropic Hookean hody to obtain the components of 
strain. 

e ■ e ■ e ■ y ■ 0 (2.4) 

x y z 'xy 

Therefore, the strain -displacement formulas reduce to 


3u 3v _ n 

H 3 ? 3 ?" 0 

and 


(2.5) 


u - f(x,y), v * g(y,z), w - h(x,z) 


where f, g , and h are continuous functions yet to be determined . 
From this observation, we conclude that cross sections do not 
distort in their own planes. In other words, the angle between any 
two lines on a cross section is not changed during the deformation 
of the beam. This makes it possible for us to define the in-plane 
displacements of any point on the cross section in terms of the 
angle of twist of a of a straight line on the section drawn from 



the x-axis to the point, as shown in Fig, 2.1. Examining the 
geometry of Fig. 2.1, we find 


u * - 0zy and v = 0zx. 


( 2 . 6 ) 


Referring to equations for stress-strain relations, we have 


9t 


ay. 


xz _ p xz 

G IT 


_ 3 /3w ,, 9u, n 
°17 ( 17 + 17> - °‘ 


(2.7) 


Then eq. (2.6) becomes 


3 <12.) _ y £<M = 0. 

dz z 


3x ''9z ; 


The first term in eq. (2.8) is zero, owing to eq. (,2.4) 
we have 

d 2 (0z) = Q 
dz 2 


( 2 . 8 ) 

Thus, 

(2.9) 


— = 0 = const. (2.10) 

aZ 

It follows that the twist relative to the section z=0 is 0z, and 

that eq. (2.6) may be expressed in terms of 0, x, y, and z. In 

summary, we express the components of displacement in the form 

u = -0yz 

V = 0XZ (2 .11) 

w = f(x,y) 

Thus, once 0 and f(y,z) are known, the complete displacement pattern 


can be evaluated. 
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Finally, from stress-strain relationships, we relate the shear 


in« stresses and strains to the displacements by 


xz 


OX 


xz 


G(& + -Is.) 

coc 3z 


T ■ Gy ■ G(t^ + - 1^) . 
y z ’yz 3y 3z 


yz 'yz oy 
Introducing eq. (2.11), we have 

t - G(-jp- - 0y) 

XZ Ox J 

X - G (-1^ + 0x) . 

yz oy 


( 2 . 12 ) 


(2.13) 


(2.14) 

(2.15) 


Differentiating the first of these equations with respect to y 
and the second with respect to z, and substituting the result 
into eq. (2.2), we have 


a 2 a 2 

0 w , o w n 

. 2 T . 2 u 

3x 3y- 


(2.16) 


which is the governing partial differential equation for the 
warping displacement. This relationship is called Laplace's 
equation for the warping function. 

The torsion problem is now reduced to one of determining the 
four unknowns, T^, t^, w, and 0. To solve this prohlem, we have 
three relationships: eq. (2.2), which, is the equilibrium condition 

and the two kinematic conditions in eqs. (2.14) and (2,15), which 
we have written in terms of the stresses. Equation (2,16) is not 
independent, since it was obtained from eqs. (2,2), (2.14), and 
(2.15). The fourth relationship necessary to solve the problem is 


the simple static condition that T and X must result in a 
r xz yz 


twisting moment of magnitude M on each cross section. 


2.2. Prandtl's Torsion Stress Function 
In 1903, Ludwig Prandtl, the distinguished German mechanician, 
introduced a scheme whe:;eby the torsion problem could be reduced 
to one of determining a single unknown. The idea is to introduce 
a twice-diff erentiable function T(x,y), called the torsional 
stress function, which has the properties 


3r 

T _ ■ "s— 
xz ay 


and 


yz 


K 

dx * 


(2.17) 

(2.18) 


When we introduce these definitions into eq. (2.2), we obtain 

<2 - 15) 


which is satisfied by any function which is continuous through its 
second derivatives. Thus, any such, continuous function will 
automatically satisfy eq. (2.2) and, therefore, lead to shearing 
stresses which are in equilibrium. The correct solution to the 
torsion problem, however, must be a state of stress providing not 
only equilibrium but also compatible strains and displacements. 
Thus, out of the infinite number of functions T which satisfy 
eq. (2.2), we must choose those which also satisfy a condition 
of compatibility. 

To arrive at this condition, we introduce eqs. (2.17) and 
(2.18) into eqs. (2.14) and (2.15), and differentiate first yith 
respect to z and second with respect to y. Then we have 
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3 ,3I\ 3 ,3w o N 

■57^ “ G ' y) 


3 / 3T n - 3 / 3w , « \ 
U ( - " G 3j ( 37 +9, ° 


( 2 . 20 ) 

( 2 . 21 ) 


Recalling that u Is also continuously differentiable, we subtract 
the second equation from the first and find 


4+4 - -“ 9 - 

3x^ 3y Z 


( 2 . 22 ) 


This is the equation of compatibility for the problem of torsion 
of beams having non-circular cross sections. Any function T 
continuous through its second derivatives which satisfies eq. (2,221 
now automatically provides Both equilibrium and compatibility. 

Any partial differential equation of this type is also called 
Poisson's equation. 

We may visualize T as being a curved surface spread over 
the cross section of the hean>. According to the definition in 
eqs. (2.17) and (2.181, the slqpe of the surface in the y direction 
is the stress in the x direction, and its slope is in the x 
direction. In fact, if n is any direction oriented a with respect 
to the x axis, as shown in Fig. 2.1, the stress directed normal to 
n is clearly 


T cosa - T sina. 
yx xz 


The slope of f in the n-direction is, by definition. 


(2.23) 


dr m 3T_ dx _3T_ djr 

dn 3x dn 3y dn 


(2.24) 
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or,, since 


dx 

dn 


cosa and 


dn 


sina, 


— = -t cosa + T sina 
dn yz xz 


(2.25) 


Hence, the slope of the T surface in any direction is equal to minus 
the shearing stress in the perpendicular direction. 

Furthermore, since we proved that no shearing-stress components 
can act normal to the boundary of the cross-section, the slope 
of T parallel to the boundary must be zero. This is possible 
only if T is a constant along the boundary. Pictorial representation 
of Prandtl's torsion-stress function Is given in Fig, 2.2. The 
slope of f parallel to the boundary curve s is 


dT 

ds 


9T dx _3F d^ 

3x ds 3y ds 



(-n) + T m 
xz 


(2.26) 


where m and n are the direction cosines of a normal to the curve. 
Since no surface forces are present, and the right side of this 
equation is zero, F must satisfy the boundary condition 



(2.27) 


Hence, F must have a constant height along the boundary. Therefore, 
without loss in generality, we assume that T is zero everywhere 
along the boundary of the cross-section. 

Any solution to eq . (2.22) provides both equilibrium and 
compatibility for cross sections of any shape. To ensure that 
it also leads to stresses which satisfy static boundary conditions 
at the ends of the beam, we must also relate F to the twisting moment 
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developed on each section. The moment developed by the shearing 
stresses about the z-axis must be 


M = M * 
z t 


(-T y + T x)dxdy 
yz^ viz J 


( 2 . 28 ) 


where the integration is carried out over the entire area of the 
cross section. Introducing P, this equation becomes 


M t * -jf ( i7 y + ir x)dxdy - 


( 2 . 29 ) 


We can now write the integral In the form 


-1 


l- B 

f d£ 


ydy 


dx - 


D 

f dr 

.dx 


xdx 


dy 


( 2 . 30 ) 


where the limits A and B stand for boundary points along a line 
y - constant, where C and D stand for boundary points along some 
line z = constant. Integrating by parts, we find 


dy. 
( 2 . 31 ) 

Now T., r„, r_, and P denote values of P at the boundary points 
A B C D 

which, according to our earlier discussion, are zero. Hence, 
we have 



1“ m 

* 


f 


\ - -f 

1 

5* 

s-/ 

1 

Pdy 

dx - 

\ 

(r D X D - I ’C X C ) - J rdx 


M 


t 


2 


P dxdy . 


( 2 . 32 ) 


A 
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This final result states that the total twisting moment on any 
section is equal to twice the volume under the surface T. 

Now we can obtain the twisting moment in terms of three 
independent quantities; G, the modulus of rigidity, which depends 
upon the material; 6, the angle of twist per unit length; and 
the constant J, which depends upon the geometry of the cross section. 


M * GJ0 


(2.33) 


J is called the torsional constant of the beam. The product GJ 
is called the torsional stiffness of the beam. The formula for 
J follows directly from eq. (2.32). 


J 


_ 2 _ 

G0 


j jrdacdy 

A 


(2.34) 


Now we can introduce the modified torsion-stress function, 


T 


xz 



(2.35) 


yz 


-Gfl il 
_Ge 9x 


(2.36) 


Performing the same algebraic operations we have done in derivation 
of Prandtl ! s stress function, we have 


2 ? 

, 81 _ 0 

n 2 ^ 2 ' 

3x 3y 


(2.37) 


This is the modified equation of compatibility for the problem 
of torsion of beams having Tion-circular cross sections.. Sometimes 
it is more convenient to use the above equation in practical 
application. 



The development of the torsion equations via stress and 
strain approaches is summarized in Table 2.1 


2.3. Solution of the Torsional Problem of the Beam Having 
an Arbitrary Cross Section by Energy Methods 
We have seen in section 2.1 that the solution of torsional 
problems is reduced in each particular case to the determination 
of the stress function satisfying the differential equation (2.22) 
and boundary conditions. In deriving an approximate solution of 
the problem, it is useful, instead of working with the differential 
equation, to determine the stress function from the minimum 
condition of a certain integral, which can be obtained, from 
consideration of the strain energy of the twisted bar. The stress 
function must satisfy the differential equation 

£r + ifr 2oe (2.38) 

3x Z 3y 

where the boundary condition is 


3r dy , _3T_ dx df _ 
3y ds dx ds ds 

For the strain energy of the twisted bar 


0. (2.39) 

per unit length, we have 


U - — 


1_ 
26 J 


(T 2 +t 2 ) dxdy 

xzyz 


HI 


<211)2 +( 3T)2 

Sx' { ‘dy J 


dxdy 


(2.40) 


If we give to the stress function T any small variation 5I 1 , which 
vanishes at the boundary, the variation of the strain energy is 
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OF POOR QiMLfTY 


Table 2.1. Development of the torsion equation, 


■ 


Two strain-displacement equations 


Two Hooke’s law equations 


One equilibrium equation 


Five equations 


Y - I* - 9y 
'xz 3x * 

3w 

Y ■ —=• + Ox 

/ z 3y 

T "Gv 
xz 'xz 

T - GY 
yz yz 

9t 3t 

+ ~3jT “ 0 

T , T , Y , Y , w - 
xz’ yz’ xz* 'yz’ 

all functions of x,y only 


One compatibility equation in 
terms of stress 

3t 3t 

-if - • - 2oe 

One equilibrium equation 

3t 3t 

xz . yz _ n 
__ — + y — “ 0 
3x 3y 


T , T 
xz’ yz 


Two stress-displacement 
equations 

T y2 ■ g< 17 + **> 

One equilibrium equation 

3t 3t 

+ m 0 

T , T ,W 
xz’ yz* 


- -2G0 


V w - 0 


3r sr 

where -5— ■ -X , -5— 

3x yz* dy 


B . C . T-0 on the boundary 


- y0 Vx + + xejpy - 0 


on the boundary 

Where Vx f Uy are the direction 
cosines of the outward normal 
of the boundary 


2 Tdxdy 


K * G l I |x“9+y"0+x ^7 - y|f|dxdy 
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J L_ 
2G 


( § 2+ ( f 2 


dxdy 


and the variation of the torque is 


‘II 


2 6Fdxdy, 


Substituting the above equations into eq. (2.40), we get 


i-JJ 


3F ^ 2 


dxdy = 20 


SFdxdy 


¥ 


€> 2 * <f > 2 


-2G0r}dxdy . 


The true expression for the stress function F is that which makes 
the variation of the integral, zero. 


r<i 

x 2 


(1L) 2 + (1L) 2 

W + l 9y' 


-2G0r }dxdy 


(2.41) 


In the approximate solution of torsional problems, we replace 
the above problem of variational calculus hy a simple problem of 
finding a minimum of a function. We take the stress function in 
the form of a series, 


r - Vo + d i r i + d 2 r 2 + 


d r 

n n 


(2.42) 


in which F , F, , r_, ... , F are functions satisfying the boundary 
o’ 1’ 2’ n 

condition. These functions are vanishing at the boundary. In 
choosing these functions, we should be guided by the membrane 
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analogy and take them In a form suitable for representing the 

function T. The quantities d Q , dp •••» are numerical 

factors to be determined from the minimum condition of the integral 

(2.41). Substituting the series (2.42) in this integral, we 

obtain, after Integration, a function of the second degree in 

d > d, , d„, .... d . The minimum condition of this function is 
o ? l 2 n 


•3U 

o 


0 


3U 

"5d^ 



(2.43) 


In this way, we obtain a system of linear equations from which, 
the coefficients d Q , dp d£, ..., d^ can be determined. Ey 
increasing the number of terms in the series (2.42), we increase 
the accuracy of our approximate solution, and by using an infinite 
series, we may arrive at an exact solution of the torsional 
problem. 

2.4. Energy-Method Solution of Torsion of the Eeams. 

Having an Airfoil Cross Section 

The cross section of conventional airfoils used in low-speed 
aerodynamic design can be approximated as follows (see Fig, 2.3). 
The upper curve can be determined by 

y “ dK(|) 

and the lower curve can be determined by 

y - -V<f> 

<(.-) « K(th) - (t'h)“[l-(th) p ] q * 

c * 


where 


(2.43) 
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The houndary conditions will be satisfied if we take for the stress 
function as approximate expression, 


f • A(y-dK) (y+d^K) . (2.44) 

Substituting equation (2.44) into the Integral (2.41) and minimizing 
the strain energy 

£-0. (2.45) 

We define the constant A 


A * - 


G0 


l-ta(d 2 +d 1 2 +dd 1 )/c 2 • 


where 


K 3 (dic/d(th)) 2 d(th) 


d-i 


t 


3 

< dt 


Using the equation already derived for the torque 

M t - 2jjrdxdy 


we obtain 


c(d+d.) 


M t • " A — T 


,3 1 

“I" 


d(th) . 


(2.46) 


(2.47) 


(2.48) 


In the case of a symmetrical airfoil, the upper and lower curves 
are described by the same expression, and the following assumptions 
are taken: 


a * 2 ’ P * 3 * 1 » d * d ]. 


± d<(f) - ±d/i tl-(|)l. 


y 


(2.49) 
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To define the constant d, we need the maximum thickness of the 
airloil 2b. First, we are going to define the position of maximum 
thickness of the airfoil from eq. (2.49) 

, .d/T-d/cfj 3 ' 

* c c 


dy _ d 

* ?^?7* 2c 


4. H 3 1/2 

+ d x 


(2.50) 


From the condition 


&•« 


c 2d 


which is the position of maximum thickness, the maximum thickness 
with respect to b is calculated using eq. (2.49). 


y ■ + d 


/f I-<|> 


to 


b . d 


d - , 


Now, we can define constant A as 
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and finally 


M t » 0.0736 


Gcd* 


i + ±MI) 

13 v 2 ' 
c 


0. 


(2.52) 


2.2 

Scaling both sides of equation (2.52) with factor — - ■ . and 

ira b 

Introducing the substitution c ■ 2a, we have the scaled torsional 
stiffness per unit length: 


K* - 0.0736 

L 1 i j* 


1 + ^(- -■-) 
13'. 2 ' 
4a 


2.2 

a +b 

„ 3,3 

Tfa b 


Then, the moment-angle of twist for a beam of an arbitrary length 
can be found as, 


M 


t 


* tta 3 b 

C *Ca 2 +b 2 


0* 


(2.53) 


2*5. Introduction to Weighted Residual Methods 
The solution of Boundary-value problems can very often be 
achieved by forming a corresponding variational problem* Under 
such circumstances , variational calculus methods are very effective 
for obtaining an approximate solution, A variational principle 
is an integral expression (.a functional) that yields the governing 
differential equations and nonessential boundary conditions of 
a problem when given the standard treatment of the calculus of 
variations. In areas other than mechanics of solids, it is more 
likely that a variational principle may not be known or may not exist. 

The governing differential equation and nonessential boundary 
condition of the problem are 



Bu * g in region V 
Lu * f on boundary S . 


(2.54) 


The exact solution u * u(x,y) is unknown. We seek an approximate 
solution, u. It may be a polynomial that satisfies the essential 
boundary conditions and contains undetermined coefficients 
a l* a 2* ••• * a n* ^us, u " u(.a,x) . We must find the values a^ 
such that u and u are ,f close M in some sense. If u is substituted 
into eq. (2.54), we obtain residuals R^ and R^ because u is not 
exact. Residuals are functions of both.* and the a^. 


^ - F^Ca.x.yl * BuOc,y). ~ gGc,yl 
m R B Ca,x,yl - Lu(x,y) - JfOc.y) 


(2.55) 


where R^ is an interior residual and 
Residuals vanish only for the exact solution, u(x,y) m u(x,y). 

We presume that u(;x,y) is a good approximation if the residuals 
are made small. This can be done by various schemes known as 
the collocation, least squares, least-squares collocation, Galerkin, 
and subdomain methods. 


R_ is boundary residual. 


2 .6. The Galerkin Method 

This method provides approximate solutions to differential 
equations directly and is applicable whether the transformation 
into a variational problem is possible or not. Having wider scope 
than variational calculus methods, it is more attractive in practice 
since there is no need to evaluate the functional even in those 
situations where it existed. 


Galerkin's method can be described as follows. Suppose we 
want to solve a linear partial-differential equation 

LIu(x,y)J = f (x,y) (2.56) 

in a region over which f(x,y) is prescribed and the boundary 
conditions are linear and homogeneous. L stands for a linear 
differential operator. 

We can rewrite equation (2.56) as 

L(u(x,y)J - f(x,y) = 0. (2.57) 

If, further, the approximate solution u(x,y) were expressihle in 
the form of a complete series of functions 

00 

u(x,y) = Z a i F 1 (x,y) (2.58) 

i=l 

satisfying the required boundary conditions, then the "exactness" 
of the solution could be expressed by the statement that in the 
region, the left-hand side of eq. (2.57) is orthogonal to every 
term in the series above. That is. 


j j (L[u(x,y) ]-f (x,y) }F i (x,y)dxdy =0, i *» 1,2,3, . . .n. (.2,59) 
S 

Since we truncate the series, eq. (2.58), to a finite number of 
terms, n, then we are in a position to use the above ideas to 
impose n-conditions of orthogonality. Thus, 


r 

s 


n 

2 a F (x,y) 
i=l 1 x 


-f (x,y) >F k (x,y)dxdy = 0 


( 2 . 60 ) 
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with i - 1, 2, 3, . n, provides a kind of averaging basis for 
evaluating the n unknown a^’s such that the approximate solution is 

__ n 

u (x,y) - 2 a.F. (x,y). (2.61) 

n i-l 1 1 

We select the "weight functions" W^ and set the weighted averages 
of the residual to zero; for i ■ 1, 2, 3, . n. Applying this 
to eq. (2.55), we have 

R i " |w i Cx,y)R L (a > x,yldV + ^(x.yjRgCa.x.yidS * 0. (2.62) 

V S 


Weight functions W^ are, by definition, coefficients of generalized 

9i7 

coordinates. Thus. W. ■ ^ , In structural mechanics, the residuals 

i o a^ 

are proportional to forces or moments, and the can be regarded 
as virtual displacement or rotation. Each integral in eq, (.2,621 
represents virtual work, which should vanish at an equilibrium 
configuration. Usually, in the problems of linear theory of 
elasticity, the Ritz coefficients are identical to the cqefficients 
found by the Galerkin method for the same system of coordinate 
functions 


2.7. Application of the Galerkin Method to the Torsion of 
Beams Having an Airfoil Cross Section 
In section 2.2, we showed that the torsion problem for a beam 
having an arbitrary cross section can be reduced to the problem 
of finding a modified torsional stress function ¥ which vanishes 
on the boundary and which satisfies the differential equation 


V 2 '? + 2 - 0 


(2.63) 


1 
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everywhere within the boundary. When this function has been found, 
the components of the shearing stress are given by 


xz 


3y 


x - -oef 

xz 3x 


(2.64) 


where 0 is the twist of unit length (in radians) and G is the 
modulus of torsional stiffness of the homogenous and isotropic 
material. Also, the torsional stiffness of a beam having a unit 
length is given by 


- 2G jj’i'dxdy 


(.2.65) 


where the integral extends over the whole section. The function ¥ 
is proportional to the deflection of the membrane in the familiar 
membrane analogy of torsion. 

The Galerkin method provides a means for the approximate 
solution of the differential equation (2.63) with the assigned 
boundary condition. Let F^, F£, F^, ... , F^ be a sequence of 
linearly independent functions of x and y which all vanish, on 
the boundary, and put as an approximation 


* “ \ a i F r 


(.2,661 


Let the result of substituting this expression on the left-hand 
side of eq. (2.63) be R^. Then is the error or residual 
in the differential equation corresponding to the chosen function 
Hf. # and it is a linear function of the coefficients a^. In 
accordance with the Galerkin method, these coefficients are 


i 
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determined by the condition that the n equations typified by 

| RgF^dxdy = Q (2.67) 

must be satisfied. This equation can easily be interpreted in 
relation to the membrane analogue. If V is proportional to the 
displacement of the membrane in a direction perpendicular to its 
undisturbed position, then is proportional to the external 
load per unit area which is left unbalanced by the tensions in 
the membrane, and equation (2,671 expresses the vanishing of the 
work done by the unbalanced loads in a virtual displacement of 
the membrane proportional to F^. Thus, the coefficients a may 
be regarded as generalized coordinates for the membrane, and the 
Langrangian equations of equilibrium are then typified by eq. (2.67}, 
Consideration of the choice of the function remains. When 
the section is symmetrical about axis x and has a smooth, boundary 
given by 

y = + th (.2*68}, 

where th.is a known function of x, the following functions, with 
p and q as positive integers, will be suitable; 

F = (th^ - y^)x^. (2.69) 

. ' pq 

These functions vanish all over the boundary. The same functions 
will serve if the boundary has a sharp corner or corners, except 
when part of the boundary is perpendicular to axis x. 

It must be pointed out that an exact solution could, in 
general, only be reached by the employment of a doubly infinite 
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i 
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set of functions However, it is found that in most cases an 

excellent approximation is obtained by the use of only two or 
three functions, and as a rule it is sufficient to use only the 
functions which are of the second degree in y, at least for narrow 
sections. 

The cross section of symmetrical airfoils can be approximated 
by the cubic oval given by 

th 2 - g 2 (2a)x(l - |) 2 (2.70) 


where c ■ 2a is the chord and g is the thickness parameter equal 
to 2b/c. This boundary very closely resembles some of the 
symmetrical airfoil sections in current use for conventional 
propeller and rotor blades. The maximum thickness occurs at one 
third of the chord from the nose, and its magnitude is 8ag/ v277 
This particular case has been worked out in detail by the thickness 
parameter method. The expression for the stiffness ohtained in 
this way is 


K 

t 


256G(2a)V 

3465 


11 2 
13 



379 

221 



(2.71) 


correct to the 6th power of g within the bracket. In applying 
The Galerkin method, the following approximation of the modified 
stress function will be adopted: 


4» . (t 2 - y 2 ) 


C 1 +C 2& 


(2.72) 


Hence, 


- V 2 f + 2 • E o + E l <?> + E 2<!>- 2 


(2.73) 
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where 


E q - 2 - c 1 (.2 + 4g 2 ) + 2c 2 g 2 

(2.74) 

E x - 6c lg 2 - c 2 (2 + 12g 2 ) 

C2.75) 

and 


h “ 12 G 28 2 - 

(2.76) 

The two Galerkin equations are 


||R a (th 2 - y 2 )dxdy - 0 

(2.77). 

S 


H^Cth 2 - y 2 )_(— Ldxdy ■ 0. 

(2.78) 

After substitution and reduction, these becjme 


Cl U3 + llg 2 ). + c 2 (5 + 3g 2 ) - 13 


^(.253 + 153g 2 ). + c 2 (li g + 81g 2 ). - 255. 

(2.79.x 

The last equations yield 


c L A- 17 + 18g 2 


c 2 A- 51g 2 

(.2.80) 

where 


A - 17 + 52g 2 f 27g 4 . 

(2.81) 


Hence, by equation (2.71)„, the expression for the torsional 


stiffness is 
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K 256G(2a)V 
t “ 3463 


221 + 489g 


221 + 676g 2 + 351g 4 


(2.82) 


The fraction inside the bracket can be converted into an ascending 
power series In g, when g is small. The erpans.iqn for the term 
in g^ can be expressed as 


11 2 
13 


1 ■ g” + g ' - oTT 8 


4 379 6 

221 


(2.83) 


The torsional stiffness of the cylinder of unit length, is. 


K _ 256G(.2a)V 
t 3465 


11 2 . 4 379 6 


1 - II 8 + 8 “ 


221 


8 


(.2.84) 


a 2 +b 2 

If we scale the expression in eq. (2.84 1 with the factor — j-j , we 

ira b 

have modified torsional stiffness In non-dimensional form 


v * _ 256G(.2a) 4 g 3 
t " 3465 


11 2 . 4 379 6 


1 “13 8 + 8 “ 


221 


(2.85) 


If the simpler approximation 


Y - ^(th 2 - y 2 ) 


( 2 . 86 ) 


were employed, then c would be found from eq. (2.80) by omission 
of the term in c 


C 1 “ 2 (2.87) 

13 + llg* 

In a similar manner, the torsional stiffness of the beam of unit 
length, is derived by 


* _ 256G(2a)V 
K t 3465" 


l— 11 

L 5 

13 + llg" 


( 2 . 88 ) 
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a 2 b 2 

Scaling the expression above with the factor - j , we get 

TTa b J 


256G(2a) 4 g 3 

3465 


(—12 

' 2 
13 + llg 


2^.2 

a -fb 

3.3 

TTa b 


(2.89) 


Finally, the relation torque-angle of twist for a beam of ur»it 
length can be established as 


M t - K t 0. 


C2.90) 


This expression can he generalized easily for any cylinder of 
length l and for a given chord 2a of the cross section: 


3 3 

M ■ if* ®. . b 

m t K t a 2^.2 

a +h 


(2.91) 


2,8. Application of the Galerkin Method to the Torsion of 
Beams Whose Cross Section is Bounded by Parabolic Curves 
A thin blade has a symmetrical section bounded by two paraholic 
(see Fig. 2.4) giving the thickness at disctance x from the center 
by 

th. m 2b [1 - (2x/2a) 2 ]. (2.92) 


The chosen function 

r - a 

satisfies the torsion boundary condition V * 0 on the boundary. 
Using the approach given in the previous section, the torsional 
stiffness per unit length is found to be 


2 2 2 
2 .th „ xV 

y - -r* 1 - 

a 


(2.93) 


K 

t 


0.1524(a) th 3 

3 3 

1 + 1 . 333 (■$—)' 


(2.94) 
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a 2 +b 2 

Scaling eq, (2.94). with the factor — s-r , we have the modified 
torsional stiffness in non-dimensional form. 


K 


0,1524 (2a) th 3 a 2 +b. 2 

~|ira 3 b 3 


(2.95) 


1 + 1.333(||) : 


Having given the geometry of the cross section and knowing the angle 


of twist per unit length, we can expand the expression above for 
a cylinder of an arbitrary length &. 




if* 0 

t l 


_ 3. 3 
Tfa b 

2 

a +b z 


( 2 . 96 ) 
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Chapter 3 

EXPERIMENTAL METHODS AND TECHNIQUES FOR THE DETERMINATION 
OF THE TORSION OF BEAMS HAVING AIRFOIL CROSS SECTIONS 

3.1, Introduction 

There are several analytical methods for the determination 
of the torsional stiffness of a cylinder having an airfoil cross 
section. Some of these methods were given in detail in the previous 
chapter for the types of cross sections which are reasonably good 
approximations of a conventional symmetric airfoil section. The 
contour lines of the airfoil section are usually presented as 
curves described with fairly complicated expressions. In trying 
to find the exact solution for a particular section by any method 
described in Chapter 2, one faces serious problems in choosing 
the appropriate stress functions and in matching given boundaries. 

In this study the experimental method is applied to three 
aluminum models with three different cross sections. Since the 
cross sections of all three models were chosen from the family 
of NACA four-digit wing sections, the approximation taken in 
analytical methods is very close to the real cross section. The 
experimental method is universal and can be used for any type of • 
beams having the airfoil cross section restricted only to the 
linear elasticity cases. 

3.2. Method of Experimental Measurement 
Since we do not know the exact position of the shear center 
of the models having airfoil cross sections, the relation, between 
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torque and angle-of-twist cannot be defined in a conventional way# 
During modeling of this system, the following assumptions were made: 

1. For models made of isotropic materials, in the elastic 
region (small deformations) t the applied torque, the angle of twist, 
the force, and the linear displacement are related linearly. 

2. The continuous beam (model) with applied torque at 
the tip is modeled as single-spring-mass system. 

3. Since the models are of small size, the weight of the 
model is neglected. 

The model was clamped to a vertical beam of much larger 
mass and stiffness than that of the model. A rigid light bar 
with, two positions for applying load was firmly attached at thn 
tip of the model. By applying the same weight on two different 
moment arms successively, the m<M*l was deflected at two different 
angles of twist. The first observation was made in the deformed 
state of the model, when the weight was applied at position 1. 

From the theory of vibrations, we know how to model some 
distributed systems as single-degree-of -freedom systems, thus 
reducing continuously distributed parameters to equivalent single 
parameters lumped at one point. Using the elementary relations 
for work done by force and torque in linear solid mechanics, we 
can model a cantilever weightless beam of length l and with given 
physical properties El as a single-spring-mass system where; 

2 

1. E is the Young modulus of elasticity given in [lb/in ]. 

2. I is the sectional moment of inertia around the axis 

4 

of bending given in in . 



A very useful relation can be established using the force-displacement 
relation for the linear spring P ■ K^A, where K^Ilb/ InJ is the 
linear spring stiffness (constant). Observing Fig. 3.1, we have 


P 



A ■ 1 




(3.1) 


In the equation above, is the equivalent spring stiffness 
of a continuous weightless cantilever beam of length. I and given 
El. . 

The cantilever beam subjected to the force and moment in 
the plane of the tip cross sectlqn is modeled as the single-spring- 
mass system. The geometry and load conditions are given in Fig. 3.2. 
Force of intensity P is applied at point 1 at the distance JL from 
shear center in the plane of the tip cross section. Point 2 
is at the distance ^ from the shear center. The difference 
l -SL is denoted as d . In Figure 3.2(,a) the model is shown in 
the displaced state, and in Figure 3.2(b) the single-spring-mass 
model is shown in the displaced state when force P is applied 
at point 1. The displacements are measured at points 1 and 2. 
Observing Figure 3.2, we can say that 

A is the linear displacement of the tip of the cantilever 
beam subjected to the force P a,t the tip. 

A^ is the linear displacement at point 1 due to the force 
P at point 1. 
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Figure 3.2. Cantilever beam model in displaced state when force P is 
applied at point 1: (a) cantilever beam model; 

(b) equivalent single-mass-spring model. 
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A^ is the linear displacement at point 2 due to the force 
P at point 1. 

Using the geometry from Figure 3.2 and equations given in 
this chapter, we can establish, the following relations; 



(3.2) 


Mi-P 


,5, a 

*1 


K_0, 

II 


(3.3) 


or 


*1 ’ 


Again, using the geometry in Figure 3.2, the following relations 
can he established; 


A 

A 

A 

A 


11 

11 

21 

21 


A + ^ 
Ml 

A + *1 Kj 


a + i z e 1 

«i 

a + h ^ • 


(3.4) 


(3.5) 


Now the cantilever beam model is observed In the second 
position. Force of intensity P is applied at point 2. Again, 
the geometry of the cantilever beam model and its single-^spring- 
mass system analogy is given in Fig. 3.3. Using the geometry in 
Fig. 3.3, the following relations are established: 


or 


M 2 - p-* 2 


K T®2- 


S' 


(3.6) 
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Figure 3,3. Cantilever beam model *in displaced state when force P is 
applied at point 2: (a) cantilever beam model; 

(b) equivalent single-mass-spring model. 
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a + Mo 
2 2 

M 2 

* + *1«2 


Now we can establish the following relations: 

P* 2 PJI A, PA, 

*22 -*21- i+ -KT- W + -^> -TT* 


A 22 " A 21 “ “k^T d a* 


In a similar manner, 


PM, Po, 

*12 '*11 ' * + - iq -- ( -* + - k 7 > 


A 12 - A 11 * d a- 

Suhstracting eq. (3.10) from eq. (3.9), we have 


A 22 " A 21 " (A 12 " A ll i " (A 2 “ V 

and f inally . 


T A 22' A 21- A L2 +A L1 


which is the torsional stiffness of the model of the length a, 

In order to get the torsional stiffness per unit length, we have 
to multiply the expression above by a: * 

Pd/ a 

t A 22" A 21“ A L2 +A L1 


(3.12) 
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If the chord of an airfoil cross section is given as 2a, we can 

2 2 

scale the expression above with the factor a -- _ 

ffa h 


C [A 22" A 2l" A 12 +A ll J 


2 ,.2 
a +b 

_ 3. 3 
TTa b 


(3.13) 


where 


* 2,2 

K m K — 

K t K t _ 3. 3 
iTa b 


(3.14) 


We can now find the torsional stiffness for any b,eam of given 
length £ and chord 2a for the sene family of airfoil section: 


K i 

T 


* * 3.3 

ira b 

K t • 


Ke 2 +b 2 ) 


(3.15) 


In a similar manner, the relation of the torque-angle of twist 
is derived; 




* TT 3. 3 

r ffa b 


t or 2 ^.U 2 \ 

Jl(A +b ) 


0 . 


(3.16) 


3.3. Experimental Technique 

The experimental technique is based on the approach presented 
in the previous section. To obtain satisfactory accuracy in 
measurement, the system for measurement is set up as shown by 
Fig. 3.4. 

The model is clamped to a vertical solid beam of much larger 
mass and stiffness. At the tip of the model, two light aluminum 
rods are attached. One of them has two attaching points for 
weight. These are shown in Fig. 3.2 as points 1 and 2, and are 
located at the distance d^ * 10{inJ, In order to avoid 
error in measurement, another bar is used as a reference line when 


I 
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Figure 3.4, 


Side view of the measurement system used. 
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measuring che displacements at points 1 and 2. Dial gages 
are used for displacement measurement, with a scale of Q..Q.01 [inj 
per division. Since the shear modulus of aluminum is relatively 
low, a small angle of twist is used on the model in order to avoid 
plastic deformations. Measurements are conducted in the 
following way: 

1, The weight is hung at position 1 and the displacements 

are measured at position 1, A^ , and at position 2, A^. 

2. The same weight is hung at position 2 and the displacements 

are measured at position 2, A ^ » and at position 1, A^, • 

The same procedure with, different weights is repeated and the 
displacements are measured. The same procedure is then repeated 
for all three models, selecting the appropriate amount of weight 
to be used in order to avoid plastic deformation. 
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Chapter 4 

PREDICTED AND EXPERIMENTAL RESULTS OF TORSIONAL 
STIFFNESS OF BEAMS HAVING AIRFOIL CROSS SECTIONS 

4.1. Introduction 

In this chapter, the predicted results are given for a group 
of models having cross sections similar to that of the symmetrical 
airfoil. Experimental results are obtained for the three models 
having NACA four-digit symmetrical sections. The results are 
given In Tables 4.2, 4.3, and 4.4. 

4.2. Comparison and Discussion of Predicted and 
Experimental Results 

In Chapter 2, the formulae are derived for torsional stiffness 
* 

per unit length for beams having cross sections close to that 
of the symmetrical airfoil section. The models have widths and 
thicknesses corresponding to the dimensions of models having 
symmetrical airfoil sections used In the experimental part. The 
formulae and the results are summarized and given in Table 4.1. 

In the experimental part, three models are used. They have 
been tested, according to the procedure given In Chapter 3, in 
the structural area of the Department of Aerospace Engineering 
of The Pennsylvania State University. The models have the 
following material properties: 

1. The material of the models is aluminum. 

2. Young's modulus of elasticity is 10.5 x 10^ [lb/in^] . 

3. The shear modulus is 4.0 x 10^ [lb/in^] . 
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2 4 

4, The mass density is 2,59 [lb sec /in ], 

The geometrical a,nd the material properties are giyen in Table 4.1. 

With, the data obtained by measurement, the torsional stiffness 

per unit length and the scaled torsional stiffness per unit 
* 

length K fc are calculated and giyen in Table 4.3. 

The. results predicted for models given in Table 4,2 and the 
results experimentally obtained for the models in Table 4,3 are 

, i 

summarized in Table 4.4. Observing these results, we conclude 

that w.e have fairly good agreement for scaled torsional stiffness 
* 

K fc obtained analytically and experimentally for the same models. 
The results obtained experimentally ere slightly below the 
results predicted. This is probably due to imperfect clamping 
of the models and the way the torque is applied to the models. 

For any given model mede of aluminum qr any metal having a 
NASA four-digit symmetrical section, we can calculate the 
torsional stiffness, if w.e know the length. I and the chord 2a; 




TtaV 

i(a 2 +b 2 ) 


(4.1) 


where is predicted and measured as a function of the shear 
modulus. Generally, we can calculate the torque for a given angle 
of twist for the models aboye as 


m t - i^e. 


(.4.2) 


Since the stiffness experimentally measured is slightly below 
predicted, it is recommended that the first one be used in further 
calculation. 



Table 4.1. Geometrical and material properties of blade models used in experimental measurement 
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Table 4.2. Analytically predicted torsional stiffness of beams having different cross sections 
scaled as the beams of unit length. 
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Table 4. A Predicted and measured torsional stiffness for three models having NASA four digit 
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Chapter 5 

THE DYNAMIC MODELING OF CANTILEVER BEAMS AS UNDAMPED 
MULTI-DEGREE-OF-FREEDOM SYSTEMS 

5.1 Free Vibration 

The analysis of a structure in free vibration provides the most ' 
important dynamic 'properties of the structure, which are the natural 
frequencies and the corresponding modal shapes. The equations of 
motion are represented in a matrix form by 

[m]{x(t)} + [k]{x(t)} « {0} (5.1) 

which represents a set of n simultaneous differential equations of 
the type 

n n 

E m .x.(t) + E k x.(t) - 0, i * 1,2,... ,n (5.2) 

i»l J j«l J 

where [m] and [k] are the mass and stiffness matvices. Before 
considering properties of frequencies and modes, it is useful to 
consider some of the properties of the stiffness matrices and mass 
matrices encountered in structural dynamics problems. 

5.2. Some Properties of [k] and [m] Matrices 
The strain energy and kinetic energy can be written in the form 
of the triple-matrix product 

V = |{x> T [k]{x> (5.3) 

T 

where {x} is the n-dimensional displacement vector and (x) is its 
transpose matrix. 
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K 11 12 
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K 21 22 * 
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• • • 
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In 

2n 


. k 


nn 


(5.4) 


and 


T = y(x} T [raj {x} 


(5.5) 


where {x} is the n-dimensional velocity vector, 



m ll m 12 * 
m 21 m 22 • 


[m) 


i - 


m nl ra n2 • 


m, 


m 


In 

2n 


(5.6) 


m 


nn 


Matrices [k] and [m] for most structures are positive-definite 
matrices, that is, when arbitrary vectors {x} and {x} are chosen and 
V and T are computed from eqs. (5,2) and (5,4), the resulting values 
of V and T are positive, except for the trivial cases {x} =0 

and {x} = 0. 

In general, the kinetic energy T is always positive-definite, 
so that [m] is always positive definite. The question remains as 
to the sign properties of the potential (strain) energy and the 
associated matrix [k] . Two cases of particular interest in 
vibrations are those in which [k] is positive-definite and [k] is 
only positive-semidef inite. When [k] is positive-semidef inite, V 
can be zero (for rigid-body modes) or greater than zero (for deformable 
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modes). When [k] is positive-semidef inite, det [k] ■ 0. Since 
the det [k] vanishes, fk] is called a singular matrix. When [m] 

<r 

is positive-definite and [k] is positive-semidef inite, the system 
is referred to as a semidefinite system, and the motion is called 
undamped free vibration. Rigid-body motion is possible because 
semidefinite systems are unrestrained, that is, such systems are 
supported in a manner in which rigid-body motion can take place. 

5.3. Some Properties of Natural Frequencies 
We are interested in a special form of solution of the set 
(5.2), namely, that in which all the coordinates have the same time 
dependence and the general configuration of motion does not change 
except for the amplitude, so that the ratio between any two coordinates 
x^(t) and Xj(t), i i j , remains constant during the motion. This 
type of motion is expressed mathematically by means of a linear 
transformation that can be represented by the matrix of the system 
eigenvectors. 

(x(t)} -[$]{q(t)} (5.7) 

{x(t)> *[<J>l{q(t)} (5.8) 

Substituting equations (5.7) and (5.8) into equation (5.1), we have 

Cm] [<J>]{q(t)} + [kH<f>]{q(t)} - {0>. (5.9) 

If we assume the following type of harmonic motion for the r-th 
mode, the characteristic equation for the r-th mode can be derived 
in the following manner. 
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{q f (t) } = {A r }sin(w r t-0 r ) 

(5.10) 

{^(t)} * -{A r }a) r 2 sin(w > t-0 r ) 

(5.11) 

(-<*) 2 [m] + [k]){<j> }{A }sin(w t-0 ) - (0) 

t r l r t 

(5.12) 


Since {A }sin(w t-0) + {0} 

^ *■ <» 

(-U> r ^ [m] + [k ] ) {c|> r ) ■ (0). (5.13) 

The equation above represents the n-th order algebraic eigenvalue 
problem. For a nontrivial solution, the determinant of eq. (5.13) 
has to vanish: 

det( [k] - a) 2 [m]) - 0. (5.14) 

This is called the characteristic equation. When the determinant 

of eq. (5.14) is expanded, a polynomial equation of degree n in 
2 

U) is formed whose roots are the eigenvalues, or squared natural 

2 

frequencies, u> r . These can be ordered from the lowest to the 
highest: 

0 £ uij 2 _< u> 2 2 < . . . . < u t 2 < . . . . £ a> n 2 . (5.15) 

It is possible for several modes to have the same frequency; 
however, this condition is unusual in structural systems. Not so 
unusual is the case in which several mode frequencies are very 
nearly the same; so close together, in fact, that very great 
computational or experimental accuracy is required to distinguish 
among then. 
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5.4. Some Properties of Eigenvectors (Modal Shapes) 

A mode shape is a unique deformation shape that the structural 
system would take on if excited solely in that mode of vibration. 
Mathematically, mode shapes are represented as vectors known as* 
modal vectors. Each element of the modal vector is the deflection 
of one degree of freedom in the structure. Equation (5.13) can be 
rewritten for the r-th mode in the following way: 

[k]{<{> r ) - w. 2 {(J> r ), r ■ l,2,3,...,n. (5.16) 

The solution vector {$ } corresponding to a particular eigenvalue 
is called a mode shape or characteristic eigenvector or an eigen- 
vector, and defines a particular structural deformation pattern 
called a mode. Because eq. (5.16) is homogenous, there is not a 
unique solution, and only a ratio among the elements of a particular 
mode shape vector can be determined. The values of the elements 
themselves are arbitrary, so that if {<j> } is a solution of eq. (5.16), 
then a r {4> r ) is also a solution, where is an arbitrary constant. 

Hence, the shape of the natural modes is unique, but the amplitude 
is not. If one of the elements of the eigenvector is rendered 
unique in an absolute sense, this automatically causes an adjustment 
in the values of the remaining n->l elements, and the ratio between 
any two elements is constant. 

5.5. Scaling (Normalizing) of Eigenvectors 
Since the modal vector represents a shape rather than the 
absolute deflection of the structure, elements of the modal vector 
are normally scaled in some arbitrary manner. The most common 
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technique is to scale the modal vector so that the largest element 
equals one. This allows the deflections at the remaining degrees 
of freedom contained in the modal vector to be expressed as a 
percentage of the maximum. It is a common misconception that modal 
vectors define actual vibratory motion of a structure. The total 
deformation resulting from external excitation is dependent not 
only on the magnitude and location of the input forces, but also 
on the participation of each mode in the structure's total response. 

The process of scaling a natural mode so that each of its 
elements has a unique value is called normalization, and the 
resulting modal vectors are called normal modes. A very convenient 
normalization scheme consists of setting 

{<l> r } T [m]{<J> r } - M r , r - 1,2 ,3,. . . ,n (5.17) 

so that the is the generalized mass or modal mass in the r-th 
mode. Generally the modal mass matrix is 

[<J>J T [mH<J>] - (5.18) 

In a similar manner, the generalized stiffness or modal stiffness 
for the r-th mode is defined as 

{^ r } T [k] (4> r } - K r , r - 1,2,3,'. . . ,n. (5.19) 

The modal stiffness matrix is 

C4>] T [k] [<j>] - rui 2 Mj . (5.20) 

The following is an especially convenient normalization for 
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T 

a general system, so that the product W } (mlC^) has the units of 
unit mass. Eigenvectors obtained in a such manner are called 
orthonormalized eigenvectors. For a system having a consistent 
mass matrix (full), the normalized i-th component of the j-th modal 
vector is 

<j> 

* 1:J “ 13 . (5.21) 


For a system having a diagonal mass matrix, the i-th component 
of the j-th modal vector Is 




A 


k-1 


Vk 


2 

J 


C5.22) 


in which {<}>^}and {<f^} are any two niodal vectors and [M] is the modal 
mass matrix of the system. The normal modes may be conveniently 
arranged in the columns of a matrix known as the modal matrix of 
the system, that is. 


(J) & 

11 12 • 


$ 


In 


A A & 

21 22 * * ‘ 2n 


m 


0 ... 4 

nl n2 nn 


(5.23) 


which holds for the general case of n degrees of freedom,. 
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5.6. Orthogonality of Modal Vectors and the Expansion Theorem 
The natural modes possess a very useful property known as 
orthogonality. This Is not an ordinary orthogonality, but an 
orthogonality with respect to the Inertia matrix [m] (and also 
with respect to the stiffness matrix [k] ). The following 
calculation is proof of the orthogonality of the modal vectors 
We are going to consider two distinct solutions, cd^ (<P r ) 

2 

and w (4> }, of the eigenvalue problem presented in eq. (5.19). 
s s 

These solutions can be written In the form 

[kH<t> r } - u> r 2 [m]{4> r } (5.24) 

[k] {<j> > - 0) 2 Cm] {4> }. (5.25) 

s s s 

T 

Premultiplying both sides of eq. (5.24) by {<}> } and both sides of 

s 

T 

eq. (5.25) by , we obtain 

{$ s } T [k]{<|> r } - 2 {<)> g } T [m] {<|> r } (5.26) 

W > T [k]W } - (U 2 {<f> } T [m]{<l> } (5.27) 

r s s r s , 

If we transpose eq. (5.27) and assume that matrices [k] and [m] 
are symmetric, and subtract the result from eq. (5.26), we obtain 

(id 2 - Cd 2 ) {<t> } T [m] W } - 0. (5.28) 

r s s r 

Since, in general* the natural frequencies are distinct, (d ^ (d , 

^ S 

eq. (5.28) provides 

{<J> g } T [m] {<J> r > - 0, r + s 


( 5 . 29 ) 
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which is the statement of the orthogonality condition of the modal 
vectors. It is easy to see that the modal vectors are also 
orthogonal with respect to the stiffness matrix [k]>‘ 

} - 0, r + s. (5.30) 

♦ S i 

We stress again that the orthogonality relations shown in eqs. (5.29) 
and (5.30) are valid only if [m] and [k] matrices are symmetric. 

The modal vectors can be conveniently arranged in a square 
matrix of order n, known as the modal matrix and having the form 

ti»] (5.31) 

where [<j>] is, in fact, the transformation matrix. All n solutions 
of the eigenvalue problem can be written in a compact matrix form: 

[k] [$] - [ml [tf>] CW 2 -J. (5.32) 

The eigenvectors {$ } (r - l,2,3,...,n) form a linearly 
independent set, implying that an n-dimensional vector representing 
any possible, but otherwise arbitrary, configuration of the system 
can be constructed as a linear combination of eigenvectors. 
Physically, this implies that any motion of the system resulting 
from an arbitrary excitation can be regarded at a given time as 
a superposition of the natural modes multiplied by appropriate 
constants, where the constants are a measure of the degree of 
participation of each mode. This approach permits the transformation 
of a simultaneous set of differential equations of motion into an 
independent set, where the transformation matrix is the modal 
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matrix [<£] . The process of deriving the system response by 
transforming the equations of motion into an independent set is 
known as modal analysis. 

To prove that a set of vectors {9^} is linearly independent, 
we assume that the vectors are linearly dependent and arrive at 
a contradiction. For the set of vectors (<P r ) to be linearly dependent, 
it must satisfy an equation of the type 

n 

C 1 {*^1 ^ + c 2^2^ + **’ + c n^n^ = ^ C r^r^ * ^ (5.33) 

r*l 

where c^ (r * I,2,..,r0 are non-zero constants. Premultiplying 

eq. (5.33) by {<£ we obtain 

s 

n 

r=l c r^s^^^^r^ * (5-34) 

T 

The triple-matrix product { <p } [m]{<p } is equal to zero for r # s 

s r 

and is less than or greater than zero for r * s. Repeating the 
operation in the equation n times, for s = l,2,3,...,n, we conclude 
that eq. (5.34) can be satisfied only in the trivial case defined 
by c^ = C£ = . .. * c * 0. Then we have {$} as a linear combination 
of {f^}, {<P 2 K •••, {9> n }, with coefficients c , c 2 , c 3> .. ., c r : 

4) = c. {<f>. } + c + ... + c } i (0). (5.35) 

xx £ i, n n 

Hence, any vector belonging to the space {$} can be generated in 
the form of linear combination as shown in (5.35) . Physically, 
this means that any possible motion of the system can be described 
as a linear combination of the modal vectors. The coefficients c 

r 

are a measure of the contribution of the associated modes } to the 

r 
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motion {$}. 

c r * ^ ' (5.36) 

Equations (5,35) and (5.36) are known in vibration theory as 
the expansion theorem. This can be applied to multi-degree-of- 
freedom systems: 

]{<}>} - 0. (5.37) 

IT S 

If matrix [m] is symmetric, we can rewrite equation (5.37) as 


"AAa + m 2^2r^2s + ••• + 'AA. ■ °* c * • <5 - 38) 


for r ■ s ■ n 


*An 2 + + • ' ‘ + ' M n 


(5.39) 


or write it in a matrix notation as 


(4> } } - M . 

n n n 


(5. AO) 


We recognize as the modal mass in the n-th mode. Generally 
for the system of order n the equation above becomes 


[<f>] [m] [<j>] 


M. 






(5.41) 


which is the generalized mass matrix or modal mass matrix. Using 
eq. (5.26) under conditions r * s = n, for the n-th mode. 


(9 } T [k] {<f> } * w 2 {<f> } T [m] (4> }. 

n n n a n 


(5.42) 


] 
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In eq. (5.42), we recognize the modal stiffness in the n-th 
mode. Generally, for the system of order n, this equation becomes 




1 "l 




2 

U) M 
n n 


(5.43) 


u . 


% -*■ 


If the eigenvectors are scaled in the following way 

[$] T [m] [$] - III 

where [I] is the identity matrix, then the scaled model -stiff ness 
matrix has the following form: 

[#J T [k] [♦] -rw 2 <J- 


w. 




0) 


(5.44) 


5.7. Response of Multi-Degree-of -Freedom Systems 
to Initial Excitation 

The equations of motion are represented in a matrix form by 

[m] (x(t) } + [k] (x(t) } * {0}. (5.45) 

Equation (5.45) represents a set of simultaneous ODE's with constant 
coefficients. Using linear system theory, a general closed-form 
solution of eq. (5.45) can be shown to exist. A convenient way 
of deriving the solution is modal analysis, which requires the 
solution of a so-called eigenvalue problem for the system. 


* 

p 


i 

I 


i 


j{ 

si 


i 
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In eq, (5.45) {x(t)} is the displacement vector or the generalized 
coordinate vector x f (t) (r <= 1,2,3, ... ,n) . VJe introduce the 
following transformation matrices: 

{x ( t ) } = [4>]{q(t)} (5.46) 

(x(t)}= [<5>]{q(t)}. (5.47) 


Substituting these two equations into eq. (5.45), we have 
t*n#Hq(t)} + [kim (q(t)} = (0). 


Premultiplying the above expression by [<J>] , we have 


[<I>] T [m] [#]{q(t)> + (<J>] T [k] [$]{q(t)> = {0>. 


(5.48) 


Using the relations in eqs. (5.43a) .atkA (5.44) derived earlier, 
eq. (5.48) can be rewritten as follows: 


or 


[I] (q(t) } + r« *J(q(t)} « {0} 




rv i 


q L (t) 

q 2 ( t} 

\ 

► + 

• 

< 

q 2 Ct) 

• 



• 2 
(0 




'o' 


-< : y 


0 

V J 


For the r-th mode, eq. (5.50) yields 


(5.49) 


(5.50) 


q r (t) + U) t q r (t) = 0 


(5.51) 


where variables q r (t) are 


identified as the normal coordinates 


of the system. Using eq. (5.46), we can write the expansion 
expression for (x(t)} as follows: 
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(x(t)} - {<t> 1 }q 1 < t) + {$ 2 }q 2 (t) + ... + {\\U) (5.52) 

or 

(x(t) } - 2 {$ }q (t) * [$){q(t)>, (5.53) 

r-1 r r 

By an analogy with the free vibration, solution of an undamped 
single-degree-of -freedom system, the solution of eq. (5.49) Is 

q p (t) - C r cos(w t - 0^) , r - l,2,...,n (5.54) 

where C r and (r - 1,2,3, »..,n) are constants of integration 
representing the amplitudes and phase angles of the normal coordinates. 
Inserting eq. (5.54) into eq. (5.46), we have 

n n 

(x(t) } ■ [4>]{q(t)} - 2 )q(t) - 2 C {4> }cos(0) t - 0 ) (5.55) 

r-1 ^ r-1 

- so that the free vibration of a multi-degree-of-f reedom system 
consists of a superposition of n harmonic motions with frequencies 
equal to the system's natural frequencies, and with amplitudes 
and phase angles depending on the initial conditions. 

Letting (x(0)} - {x Q } and (x(0)} - {x Q } be the initial 
displacement and initial velocity vectors, respectively, eq. (5.55) 
leads to 

n 

{x Q } - S C {$ }cos9 r (5.56) 

r-1 r 
n 

(xq) - E c r u r {4> r }sin0 r . 


(5.57) 
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Premultiplying eq. (5.56) and (5.57) by {$ } T [ra), 

s 

{$ r } T [m]{x Q } =* C r cos0 r (5.58) 

^~{<t r } T (m]{x 0 } - C r sin0 r . (5.59) 

Introducing these two equations into eq. (5.55), we obtain the 
general expression 

^ T 

(x(t)} «E [{$ } [m]{x }costi3t 
r-1 r or 

. + {G^tm]^}^ sinei^t] {* r ) (5.60) 

which represents the response of the system to the initial dis- 

* 

placement vector (Xq) and the Initial velocity vector {xq}. Now 

we assume that the initial displacement vector {x Q } resembles a given 

normal mode, say {$ }, whereas the initial velocity vector is zero. 

s 

Introducing {xq} ■ Xq{$ 3 ) and {x^} ■ (0), the response is 

n T 

(x(t)} - E [q Q {$ K [m]{x 0 }cosu> r t]{$ r }- (5.61) 

r-1 

5.8. The Dynamic Response of Multi-Degree-of -Freedom 
Systems to a Harmonic Force 

The general response of a discrete system, i.e., the response 
to both initial and external excitations, is governed by a set of 

simultaneous (coupled) equations, but can be rendered independent 

* 

(uncoupled) by means of linear transformation that can be 
represented by the matrix of the system eigenvectors. Because 
these eigenvectors are commonly referred to as modal vectors, the 
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transformation matrix is called the modal matrix and the decoupling 
procedure itself is known as modal analysis. The independent 
ordinary differential equations resemble those of low-order systems 
and can be solved with relative ease. The important feature of 
these low-order systems is that they are described in terms of a 
single dependent variable. 

For a linear system, solutions to excitations can be obtained 
separately and then added up linearly to obtain the combined response. 
This is the essence of the so-called principle of superposition, 
a very powerful principle that applies to linear systems alone. 

In vibration, we encounter the various types of periodic excitations, 
which are not necessarily harmonic. Any periodic function, however, 
can be represented by a convergent series of harmonic functions 
whose frequencies are integral multiples of a certain fundamental 
frequency, U)q, provided that it satisfies certain conditions. 

The integral multiples of the fundamental frequency are called 
harmonics, with the fundamental frequency being the first harmonic. 
Such a series of harmonic functions is known as the Fourier Series, 
and can be written in the form 


F(t) ■ j a + £ (a cosptogt + b sinpoi^t) 
0 p*l P 15 


(5.62) 


2n 


where p are integers, p - 1,2,3,. . . ,n, and T is the period. 

The coefficients a and b are given by the formulae 
P P 


T/2 


2 

*P " T 

-T/2 


F(t)cospu) 0 tdt, p * 0,1, 2, 3, . . . ,n 


(5.63) 
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T/2 

b p = ^ j F(t)sinpUi 0 tdt, p = 0, 1,2 ,3, . . » ,n (5.64) 

-T/2 

and they represent a measure of the participation of the harmonic 
components cosp^Qt and sinpi^t, respectively, in function F(t). 

The term y represents the average value of F(t). We consider 
the equation of motion of an undamped n degree-of-freedom system 
subjected to a harmonic force, expressed in terms of Fourier 
components of the fundamental frequency only. The same procedure 
could be repeated for higher harmonics. 

The amplitudes of harmonic forces are distributed arbitrarily 
on all mass elements of the system. The exciting force can be 
written as follows: 

F 

F(t) = —■ + F os SinU) f t + F 0c CosU) f t (5.65) 


or 


F 0 mt-no 

F(t) = -f + F 0r e 


whare 


Or 


m 


+ F 

Os T r 0c 


(5.66) 


= arctan G~-) (5.67) 

£ 0c 

and where is the forcing frequency of the system and \|» is the 
phase angle. The equations of motion for an n degree-of-freedom 
system can be written in matrix form as follows: * 

j(u> t+l|>) 

[m] tx(t) } + [k] (x(t) } = (F(t)} = y {f q } + (F Qr }e 


. (5.68) 


The linear transformation matrices can be expressed as 


(x(t)> = (5.69) 

(x(t)} = [4>]{q(t)> (5.70) 


where the transformation matrix [4] was introduced earlier as the 

scaled modal matrix. Also, in a previous section we showed the 

orthogonality conditions for two modes, u> and u> , w ^ to 

s r s r 

{♦ } T [mj{# } * 0 r ^ s (5.71) 

s r 

{* g } T [k]{* r } - 0 rilsi (5,72) 

Substituting eqs. (5.71) and (5.72) into eq. (5.68), we have 

j(W f t+^) 

[m] [$] (q(t) } + [k] [#J (q(t) } * (F(t)} ■ j {F q } + {F^} e 


Premultiplying (5.73) by [4>]^, we have (5.73) 

T TIT T j (W,t-H|0 

[4>] [m] (<H (q(t) } + [4>] [kl [$] {q} — r fF 0 > +m{F 0r }e 

(5.74) 

Using relations derived in a previous section, eq, (5.47) reduces 
to 

? 1 T T j (CO t+ip) 

[I] (q(t) } + ru j {q(t)> = y[$r{F 0 > + [4>] (F^Je £ . (5.75) 


This is a so-called modal equation, which decouples the initial 
system of equations into n independent equations. The terms on the 
right side are called the components of modal force. They give 
the extent to which the components of exciting forces participate 
in every mode. 
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The equation of motion in the r-th mode due to the exciting 
forces lumped at all masses, having arbitrarily distributed 
amplitudes, evolves from eq. (5.75) as 


o It t 

q r (t) + ^ 2 q r - + {$} T {F 0r )e f 

The solution of eq. (5.76) can be written as 

q c U) ■ <q r (C)) h01n + (q r (t)) patt . 
and the homogenous solution can be expressed as 


(5.76) 


(5.77) 


^^r^hc® 


E cosU) t + F sinU) t. 
r r r r 


The particular solution can be found separately for each 
component of the forcing function in eq. (5.76): 

' ( S>p.rt 1 + ( ’r>p.rt 2 . < 5 - 78 > 

We assume the following types of solutions: 


^part 1 ' C • 


(5.79) 


Substituting eq. (5.79) into eq. (5.76) , and omitting the second 
term on the right-hand side of eq. (5.76), we have 


(q r ) part 1 
^Sr^part 2 


^ q r^part 2 


2 co 


Vv T v 


V 


j (CO f t+^) 


(5.80) 

(5.81) 


■v* 6 


2 j (u f t-H(0 


(5.82) 
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Substituting eqs. (5.81) and (5.82) into eq. (5.76), and omitting 
the first term on the right-hand side of the latter equation, we have 


^Sr^part 2 


{$r H F Q r ) j(w f t+0 


, 2 2 
to - w, 
r f 


(5.83) 


The general solution for the r-th mode can be written as 


Vt) 


E costo t + F sinto w -r — A LT 

r r r r 2u) L r 

r 


t + > T if„} + 1 V T V « J( "f t+W . 


0 J 2 2 
o> - 0J- 

r f 


(5.84) 


For convenience in further calculations, eq. (5.84) can be rewritten 
as 

q (t) - E cosW t + F sino) t + — ~ {$> } T {f 0 > + 
r rr r r 2U) r U 

r 


{VV S > {$ r }T{F 0c } 

2 T2 Sin V + ~ r~ 2 C09 V‘ 


(5.85) 


to - to. 

r f 


u» - u_ 
r f 


To define the constants E r and F , we apply the initial 
conditions given. The displacement and velocity vectors at t • 0 
are given as {x(0)}and (x(0)}. To impose the initial conditions 
from the original system to a new, transformed system, we use 
eqs. (5.69) and (5.70) at t * 0: 


(x(0) } - [41 (q(0) } 
(x(0) } - [4] (q(0) ). 


(5.86) 

(5.87) 


From eqs. (5.86) and (5.87), 
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{q(0)} * [$>] ^(x(0)} 


{q(0)} - [$]" 1 {i(0) >. 


(5.88) 

(5.89) 


In order to avoid the inversion of matrix .[<&]', we are using the 
following relation: 

[4] T tm][4] - [I]. 


-1 


Post multiplying the relation above by [$] , we have 


W 1 “ [#] T [ml. 


(5.90) 


Substituting eq. (5.90) into eqs. (5.88) and (5.89), we have 


{q(0)}» [*] ImHx(O)} 
(q(0)}- [$] T [m] (x(0) } . 


(5.91) 

(5.92) 


With time equal to zero, eq. (5.85) reduces to 

1 T [♦] T {F 0 } 

{q_(0)> - {E } + n } + — X * 

r r 2d) 0 (D L - <i)/ 

r r f 


(5.93) 


Substituting eq. (5.91) into eq. (5.93), we can calculate the 


vector of constants {E r ): 

} 

{E > - [*] T [ml{x (0)} - —• [$] T {F 0 > - —5 , 

r r 2U 2 0 to_ 2 -(i), 2 


(5.94) 


r f 


Finding the first derivative of eq. (5.85), 


V* 1 > 

tq^(t))* -{E }(d sincu t + {F }u cosia t + — r = — cos(D.t 

r rrr rrr w 2 _uj 2 f 


V*‘ {F o s > 
» 2 . Ue 2 


sino^t , 


(5.95) 
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With tine equal to zero, eq. (5.95) becomes 

W f (*] T {F n > 

(q(0) } = (F )0J + — . 

** r or - w/ 


(5.96) 


Substituting eq. (5.92) Into eq. (5.96), we can calculate the 
vector of constants {F^. } ; 

v-r t«l T hHx(0)} - V . ' . (5.97) 

r (W r -o» f ) 

Now Eq. (5.85) can be written in a final form as 


(q(t)} ■ [♦j T [m] t{x(0) JcosU^t + ~ {x(0) }sinu) r t J 

r 

1 - I$]T{F 0c } 

- — srC^l 1- {F n }coso) t + — s s- (cosOJ-t — slnu) t) 

20) 2 0 r 0) 2 - 0) 2 f r 
r r f 


1*1 {F 0s > W f 

+ — 5 2 ( cosU f t “ rr* alna)_t) , r - l,2,3,..,n. (5.98) 

“r '“f 


The solution to all n equations in eq. (5.98) must be calculated; 
then the finite displacements of masses can be obtained by super- 

l 

position of the response in each mode: 

(x(t)} - [4>] (q(t) }. (5.99) 

Therefore, in summary, the response analysis by mode super- 
position requires, first, the solution of eigenvalues and eigenvectors 
of the problem in eq. (5.32), then the solution of the decoupled 
equilibrium equations in eq. (5.98), and finally, the superposition 
of the response of each eigenvector as expressed in eq. (5.99),. 
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The essense of a mode superposition solution of a dynamic 
response is that frequently only a small fraction of the total 
number of decoupled equations need be considered in order to obtain 
a good approximate solution to the actual response of the system. 
Considering the problem of selecting the number of modes to be 
include. i in the modal superposition analysis, it should always 
be kept in mind that an approximate solution to the dynamic 
equilibrium equations is sought. Therefore, if not enough modes 
are considered, the governing equations of motion are not solved 
accurately enough. This means that equilibrium, including the inertia 
forces, is not satisfied for the approximate response calculated. 

This question arises especially when a lumped mass model is 
used to approximate a system having a uniformly distributed mass. 

There is no unique criterion to determine the number of modes to 
be taken. Usually in shock and stead v-state problems driven on 
high frequencies, "effective modal mass" is used as parameter. A 
modal effective mass can be identified in matrix form for the n-th 
mode as; 




{l} T [m]{ < J ) }{4> } T [ra]{l} 
n n 

M 

n 


( 5 . 100 ) 


or in summary form as 

(M 


effV n 2“ 


n 

E m.4> 

i-1 1 in 


X 

where ra^ is the mass of the i-th point 9 {1} is a row matrix of 

ones, and M is the modal mass in the n-th mode defined by 
* n 



Usually a good indication that enough mass has been considered by 
modal analysis is that the effective modal mass is very close to 
the actual mass of structure, that is, more than 95% of the total 
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Chapter 6 

DAMPED MULTI -DEGREE-OF-FREEDOM SYSTEMS 
IN STRUCTURAL DYNAMICS 

6.1. Introduction 

All real systems dissipate energy through one of various types 
of damping mechanisms, such as structural viscous or Coulomb damping. 
The most common method of taking into account the dissipation of 
energy in structural dynamics is to assume in the mathematical model 
the presence of damping forces which have magnitudes proportional to 
the relative velocity and directions opposite to the motion. This 
type of damping is known as viscous damping, because it is a kind 
of damping that will occur when there is motion in an ideal viscous 
fluid. The inclusion of this type of damping in the equations does 
not alter the linearity of the differential equations of motion. 

Since the amount of damping commonly presented in structural systems 
is relatively small, its effect is neglected in the calculation 
of natural frequencies and mode shapes, However, to uncouple 
the damped differential equation of motion, it is necessary to impose 
some restrictions on the values of the damping coefficients in the 
system. These restrictions are of no consequence due to the fact that 
in practice it is easier to determine or estimate modal damping 
ratios rather than absolute damping coefficients. In addition, when 
solving the equations of motion by the modal superposition method, 
only damping ratios are required. When the solution is sought by 
other methods, the absolute values of damping coefficients may be 
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calculated from modal damping ratios by any of the various methods. 

6.2. General Consideration 

The equation of motion of a linear multi-degree-of -freedom 
damped system is 

. 1 j 0» f t+« ,, n 

[m] (x(t) } + [c](x(t)> + [kJ (x(t) } - | {F Q > + {F 0r >e f 

In general, the coefficient matrices in the above equation, [m] , 

[c] , and [k] , may have nonzero coupling terms (e.g., c^ = c^ ^ 0), 
so that to solve eq. (6.1) in its present form would require simultane- 
ous solution of n equations in n unknowns. The purpose in this section 
is to outline the mode-superposition method by which such a set of 
coupled equations can be transformed into a set o£ uncoupled 
equations through use of the normal modes of the system. Equation 
(6.1) is the original set of coupled equations of motion for an n-degree- 
of-freedom system, where x(t) may be physical or generalized 

coordinates. The response of the system to the excitation force 

1 , 

■j {Fq} + {Fg^Je and to the initial conditions lx(0)} and 

{x(0)} is sought. 

The first step in a mode-superposition solution is to obtain 
the natural frequencies and modal shapes of the system. In Chapter 
5 we gave a very detailed treatment of these quantities. Since 
the structural systems can be considered as lightly damped systems, 
the most important dynamic properties, the natural frequencies and 
the corresponding modal shapes, are not practically affected by 
the presence of damping. We assume that if there are any repeated 
frequencies, the associated modes have been orthogonalized, so that 

c - 
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the orthogonality equations 

{$ r } T [m]{<t> r } - {4> r > [lei {<J> r > - {0} (6.2) 


are satisfied for all r + s. The modes are then collected to form 
the scaled modal matrix [$] already introduced in section 5.4.4, 
that is, 

[*] - [{*.}{*,}...{* >...{* }]• ( 6 . 3 ) 

X z r n 

The key step in the mode-superposition procedure is to intro- 
duce the coordinate transformation, exactly the same procedure used 
for an undamped multi-degree-of -freedom system in Chapter 5: 

{x(t)> - [*]{ q(t)> - l {* J .}q r (t) . (6.4) 


The coordinates q r (t) are usually referred to as normal or principal 
coordinates. Equation (6.4) is substituted into eq„ (6.1) and the 
resulting equation is multiplied by [$] to give the equation of 
motion in principal coordinates, namely, 


[$]{q(t)} + l*] T [c]($]{q(t)} + [$] T [k] [*]{q(t)> 


(<M 


j (to t-j-tp) 

2 {F 0 }+{F 0r> C 


(6.5) 


Using relations in eqs. (5.43a) and (5.44), eq. (6.5) is reduced to 


[I ] (q(t) } + [C ){q(t) + [K ]{q(t)} * (F (t) ) (6.6) 

S 3 3 


where 


[I] 




T 

[$] [m] [$] * scaled modal mass matrix, 

T 

[$] [c] [<t>] » scaled modal damping matrix, 


(6.7) 
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tK>] - I#] 
(f (t) } * I«] 

S 


Ik] t$] 

(I{ F } 
2 l 0 J 


scaled modal stiffness matrix, 

j(W f f,+^Tl 

{Fg^Je I - scaled modal force vector. 


Due to the orthogonality conditions of eq. (6.2), [m] and [k] 

are diagonal matrices, so the equations of motion in modal coordinates 

(eq, (6.6)) are coupled only through a nonzero off-diagonal in the 

scaled modal damping matrix [C ] . In the derivation of the modal- ' 

s 

coordinate equations of motion, it has been assumed that the normal- 
coordinate transformation serves to uncouple the damping forces in 
the same way that it uncouples the inertia and elastic forces. 

In the following sections, we are going to consider the conditions 
under which this uncoupling occurs. 


6.3.. Conditions for Uncoupling of Damping Forces 
Considering the analysis of systems in which damping effects 
cannot be neglected, we still would like to deal with decoupled 
equilibrium equations (eq. (6.6)) merely to be able to use essentially 
the same computational procedure, whether damping effects are included 
or neglected. In general, the damping matrix [cj cannot be constructed 
from element damping matrices, such as the mass and stiffness matrices, 
and its purpose is to approximate the overall energy dissipation during 
the system response. 

6.3.1, Evaluation of the Damping Matrix for Any Set of Specified 
Modal Damping Ratios 

In principle, the procedure can be explained by considering 
the complete diagonal matrix of the damping coefficients, which 
may be obiained by premultiplying and pos multiplying the damping 
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matrix by the scaled modal matrix 




Ct] T [c] [4>] 


2 ;.(!>, 0 . 

0 2C 2 w 2 




2C (0 

n n 


( 6 . 8 ) 


It is evident that the damping matrix [c] may be evaluated by 
postmultiplying eq. (6.8) by the inverse of the scaled modal matrix 
and its inverse, such that 

lc] - [$]" T [C S ] I*]" 1 . (6.9) 


Therefore, for any specified set of modal damping ratios {?}, 

matrix PC -J can be evaluated from eq. (6.3), and the full damping 
s 

matrix [c] from eq. (6.9). However, in practice, the inversion of 
the modal matrix is a large computational effort. In Chapter 5, 
wo proved that 

[♦r 1 - [$] T [m]. (6.10) 

Substituting eq. (6.10). into eq. (6.9), we obtain 

[c] * [m] [*] tC>J C$] T [m], (6.11) 


Since, in the r-th mode, C gr * 2^0)^, substituting this expression 
into eq. (6.11), we have 


[c] 


[m] 


— n f-i 

r-1 


r r 


[m] . 


( 6 . 12 ) 


The damping matrix [c] obtained in eq. (6.12) will satisfy the 
orthogonality property; therefore, the damping term in the differential 
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equation (eq. (6.1)) will be uncoupled with the same transformation 
(eq. (6.4)), which serves to uncouple the inertial and elastic 
forces. 

It is of interest to note In eq. (6.12) that the contribution 
to the damping matrix of each mode is proportional to the modal 
damping ratio; any undamped mode will contribute nothing to the 
damping matrix. As mentioned previously, the modal damping ratios 
are the most effective measures of the damping in the system when 
the modal analysis is carried out. 

6.3.2. Defining a System Damping Matrix Using Rayleigh Proportional 
Damping 

The mode superposition analysis is particularly effective if it 
can be assumed that the damping is proportional in the following 
manner : 

{* t > T Cel {* d > - [20)^3 6^ (6.13) 

where is a modal damping ratio and 6^ is the Kronecker delta 

(6y “ 1 for 1 * j, 6^ ” 0 for i t j). Using eq. (6.13), it is 
assumed that the eigenvectors {t^}, i ■ 1,2 ,3, ... , n, are a 1 so 
[cj -orthogonal. 

In considering the implications of using eq. (6.13) to take 
account of damping effects, the following observations can be 
made. Firstly, the assumption in eq. (6.13) means that the total 
damping in the structure is the sum of individual damping in each 
mode. In fact, the ability to measure values for the damping ratios 
and thus approximate in many cases and in a realistic manner 
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the damping behavior of the complete structural system, is an important 
consideration. A second observation relating to the mode superposition 
analysis is that for the numerical solution of eq. (6.6) using the 
decoupled equations, we do not calculate the damping matrix [c], 
but only the stiffness and mass matrices, [k] and [m]. 

As discussed, damping effects can readily be taken into account 
in mode superposition analysis provided that eq. (6.13) is satisfied. 
Rayleigh damping can be assumed to be of the form 

[c] - a 0 [m] + a^k) (6.14) 

where a^ and a 1 are constants to be determined from two given damping 
ratios that correspond to the two unequal frequencies of vibration. 
Applying the orthogonality' condition to eq. (6.14), that is, pre- 
multiplying both sides of this equation by the transpose of the r-th 
mode and postraultiplying by the modal matrix [4)], we obtain 

{4> r } T [c] [4>] = a 0 {4> r ) T [m] [4>] + [k] [$) . (6.15) 

The orthogonality conditions (eq. (5.30)) reduce eq. (6.15) to 

(4> r } T [c] m = a 0 {4> r } T [c]{4> r } + a 1 {4) r ) T [k]{4> r }. (6.16) 

Using eqs. (5.43a) and (5.44), eq. (6.16) is reduced to 

{4> r ) T [c] [*] = a Q + a^ 2 (6.17) 

which shows that when the damping matrix is of the form of eq. (6.14), 
the damping forces are also uncoupled with the transformation (eq. 
(6.15)) . However, it can also be shown that there are other 
matrices formed from the mass and stiffness matrices which also 
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satisfy the orthogonality condition. 

In general, a damping matrix that satisfies the relation 
presented in eq. (6.13) is obtained using the Coughey series: 


p-1 . \ 

[c] - [m3 Z a,{[mj~ [k)> 
x-0 A 


(6.18) 


where the coefficients a^, \ m 1,2,3, ... ,p, are calculated from the 
p simultaneous equations: 

"a!n 




+ a l“i + a 2 W i + ' * * + a p-l W i P 3 


i * 1 , 2 , 3 , , » » , p . 

(6.19) 


We note that with p ■ 2, eq. (6.18) is reduced to Rayleigh damping 
as presented in eq. (6.14). An Important observation is Chat if 
p > 2, the damping matrix [c] in eq. (6.18) is, in general, a full 
matrix. Since the cost of analysis is increased by a very significant 
amount if the damping matrix is not banded, in most practical analysis, 
Rayleigh damping is assumed. In practice, reasonable Rayleigh 
coefficients in the analysis of a specific structure may often be 
selected using available information on the damping from a conducted 
vibration test or on the damping characteristics of a typical similar 
structure, i.e., approximately the same ag and a^ values are used in 
the analysis of similar structures. 

By taking two terms corresponding to k ■ 0 and k * 1 in 
eq. (6,18), it is possible to compute the damping coefficients 
necessary to provide uncoupling of a system having any desired damping 
ratios in any specified number of modes. For any mode r, the modal 


ORIGINAL PAGE IS 
OF POOR QUALITY 


87 


damping is given as 

C r - {$.} fcH* r > - 2? r w r * (6.20) 

If [cj as given by eq. (6.18) is substituted in the expression for 


C , we obtain 
r’ 


C = {*_}%] 2 a.{[mJ _1 [k)} A {4> }. 

r r A=0 A r 


2 

Now, using relation [k] {$ } = to [mJ (4> } and performing several 

r r r 

algebraic operations, we can show that the damping coefficient 
associated with any mode r may be written as 

2 A 


( 6 . 21 ) 


from which 


C ■ Ea..UL * 2 C to 
r ^ A r *r r 


? = ■stt- E a, to 

^r 2to , A r 
r A 


( 6 . 22 ) 


(6.23) 


Equation (6.23) may be used to determine the constant a^ for 
any desired values of modal damping ratios corresponding to any 
specified number of modes. For example, to evaluate these constants 
specifying then modal damping ratios, C 2 > C3 » *•• ,H we may 

choose A ■= 0, 1,2 ,3, . . . ,q. Equation (6.23) gives the following 


system of equations 


V 


c; 


n 


1 

> ' 2 


L J L_ r) 


_i_ 

to 

1 

_ 1 _ 

U 2 

_1_ 

to„ 


to 


W 1 

3 

W 1 


“l 

W 2 

W 2 3 


W 2 

W 3 

“3 3 

“3 5 

W 3 






. . to. 


n+2 


to 


n+2 



r 

a i 


a 2 

< 

a 3 


lS 


r (6.24) 


88 


or, written symbolically, 

{$} - \ CQ x ]{a) (6.25) 

where [Q^J is a square matrix having different powers of the 
natural frequencies. The constants (a) could be obtained as 

{a} - 2iq x r 1 {?>. (6.26) 

Finally, the damping matrix is obtained by the substitution of eq. 
(6.26) into eq. (6.18). 

It is interesting to observe from eq. (6.23) that in the special 
ca fje when the damping matrix is proportional to the msis [c] ■ a^IM], 
(i ■ 0), the damping ratios are inversely proportional to the natural 
frequencies; thus, the higher modes of the structures will be given 
very little damping. Analogously, when the damping is proportional 
to the stiffness matrix ( [c ] * a^k]), the damping ratios are 
directly proportional to the corresponding natural frequencies, as 
can be seen in eq. (6.23) when i ■ 1. In this case, the higher modes 
of the structure will be very heavily damped. 

Rayleigh coefficients are to a large extent determined by the 
energy dissipation characteristics of the structural materials. 

In the above discussion, we assumed that the damping characteristics 
of the structure can be represented appropriately using proportional 
damping, either in a mode superposition analysis or in a direct 
integration procedure. In many analyses, the assumption of 
proportional damping (i.e., that eq. (6.13) is satisfied) is 
adequate. In analysis of a structure with widely varying material 
characteristics, nonproportional damping may need to be used. 



We should mention at this point the circumstances under which 
it will be desirable to evaluate the elements of the damping matrix, 
as in eqs. (6.21) or (6,18), If, b vs been stated that absolute 
structural damping is a rather difficult quantity to determine or 
even to estimate. However, modal damping ratios may be roughly 
estimated on the basis of past experience. This past experience 
indicates that values for the modal damping ratios in structures 
are generally in the range of 2% to 10%, probably not exceeding 20%. 
In other words, on this basis, and giving some consideration to the 
type of structure and materials utilized, we can assign numerical 
values to the modal damping ratios. The modal damping ratios are 
then used to determine the damping matrix, which is needed explicitly 
when dynamic response is obtained by some analytical procedure 
other than modal analysis, i.e., step-by-step integration of a non- 
linear system. 

6.4, The Dynamic Response of the Damped Multi-Degree-of-Freedom 
Systems Subjected to Harmonic Forces 
In Chapter 5 we gave a very detailed treatment of the dynamic 
response of an undamped multi-degree-of-freedom system subjected to 
harmonic forces. In this section we are going to consider the damped 
MDOF system. The general response of a discrete system, i.e., the 
response to both initial and external excitations, is governed by a 
set of simultaneous coupled equations. This set can be uncoupled by 
means of linear transformation that can be represented by the matrix 
of the system eigenvectors. The decoupling procedure itself is 
known as modal analysis, since the eigenvector matrix is referred to 
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as a modal matrix. We assume the proportional damping of the system; 
in this case, the classical modal matrix does uncouple the equations 
of motion. We will consider in detail the circumstances under which 
the decoupling of damping forces is possible. 

Again, the system is going to be subjected to a periodic 
forcing function, which may be expanded into a Fourier Series whose 
terms are sine and cosine functions of successive multiples of the 
fundamental frequency, as shown in eq. (5.62). The coefficients of 
a series are given by eqs. (5.63) and (5.64). The response of the 
system is then obtained as the superposition of the response for each 
term of the Fourier expansion of the exciting function. This function 
is given by eq. (.5.65). The equations of motion of an n-degree-of- 
freedom system can be written in a matrix form as: 

. i j (u) t+tp) 

[m] (x(t) } + [c] {x(t) } + [k] { (t)} - {F(t)> - y{F Q > + {F Qr }e f 

(6.27) 

Quantities {F Qr } and ip are defined by eqs. (5.66a) and (5.67). 
Introducing the linear transformation matrices using the scaled 
modal matrix ($>], we have 

(x(t)) - [<t>] (q(t) } 

(x(t)} - [$] (q(t) } (6.28) 

(x(t) } - [4>){q(t)}. 

In a previous section, we showed the orthogonality conditions 

for two modes, (0 and ul , 01 / U3 , can be written as follows: 

' s r s r’ 
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{4> g } T (ml{<I> r } * 0 t i 1 s 

{$ g } T [cJ {$.} -0 t i* s 

{$ } T [k] {♦ } - 0 r ^ s, 

S X 


(6.29) 

(6.30) 


Substituting equations (6.29) and (6.30) into eq, (6.27), we have 


WC*]{qCt)> + [cH*HqCt)> + Ik][<t]{qCt)> - 

1 j (W f t+W 

■ (F(t) } - j {F 0 } + {F Qr }e . (.6.31) 

T 

Preoultiplying eq. (6.31) by [$] , we have 

[*j{q(t)> + t*) T (c](*l{q(t)} + Wl T tk) m(q(t)} - 

, * T (jw t+ljj) 

- |t*r(F 0 } + t*I {F 0r >e f _ (6.32) 

Using the orthogonality relations derived in a previous section, 
eqs. (6,32) are reduced to 

lll{q(t)} + t2;uj(q(t)) + [U) 2 g{q(t)> - 

T rn m j ( b+l|<) 

- |[<&] {F 0 > + m (F 0r )e £ . (6.33) 


The n equations above are decoupled into n independent equations. 
The terns on the right-hand side represent the so-called scaled modal 
force. The equation of motion in the r-th mode, due to the exciting 
forces lumped at all masses having arbitrarily distributed amplitudes, 
evolves from eq, (6.33) as 


q' r (t) + 2? r U) r q t (t) + W. q r (t) * 

1 T T i (0) t+r ) 

- 2 { V <V + W r> tF 0r )e 


. (6.34) 
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The general solution of equation (6,34) can be written as 


q r (0 - <q r (t» h + (q r (t)) p . (6.35) 


The homogenous solution Is 

-Cw t 

(q r ( t )) h " e r [C r cosa) dr t + D r sln(O dr t]. (6.36) 


The particular solution Is composed of two components: 


( q r (t)) p “ <q r (t)) pl + (q r (t)) p2* <6.37) 


We assume the following type of solution: 


<q r (t)) P l “ C ' 


(6.38) 


Substituting eq. (6.38) Into eq. (6,34) and omitting the second term 
on the right-hand side of the latter equation, we have 


^4 » r > V- 


( «r (t,> p2 ’ It 1 


i (w f t+^ r ) 


j (u f t+ip r ) 


( q r < t )) p 2 " J“ f Qe 

2 i <01 t+T|i ) 

( \ (t)) p2 ■ -“f <> e 


(6.39) 


(6.40) 


(6.41) 


(6.42) 


Substituting eqs. (6.40), (6.41), and (6.42) into eq. (6.43) and 
omitting the first term on the right-hand side of the latter 
equation we have, 
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< q r< C » P 2 


($ r } {F Qr } j (u f t -hp r ) 


W r - w f + J 2? r ? f W r 


(6.43) 


Rewriting the complex number in equation (6.43) in terms of the 

magnitude and phase angle, eq. (6.43) becomes 

T i t+i|» ) 

{F 0r>* f r 

{q r (t)} p2 ‘ -y======== - — , 

2 )+(2c (o f U) ) 2 e j 


or 


t< >r< c) >p2 


{* r ) T (F 0r }e * r 1 


/ ( u r 2 - 0 ) f 2 ) 2 +( 2 ?i . M£mr > 2 


(6.44) 


where 


6 “ arctan 

r 


2 C r w f w r 


2 2 
U) -air 

r f 


(6.45) 


Let ^ -<S 
r r 


® . Then eq. (6.44) is reduced to 




p2 


<r> i (u) F t:-hp ) 

<V g 0 ,> e ‘ 

/(“ r 2 -“ f 2 ) 2 +(2? r u f u r ) : 


The general solution for the r-th mode can be written as 


q r (t) 


•C lil t 

*r r 


CC r COS “dr t * D r sina) dr tJ + ~^T ^r }T{F 0 } + 

2u)_ 


T J(u f t+0 ) 

, tyVo,}. 

v /w r 2 -t0f 2 )+( 2 ?r“fW r ^ 2 


(.6.46) 
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The equation above can be rewritten for a very general case as 


n -5 w t n it. 

{q r (t) } « E e r r [C^coslO^t + D^inw^t] + E — ^{^KCFq} + 


T J CU). t+0 ) 
n (* } X (F n _}e f 
+ Z r ur 


'- 1 


(0) -0). ) +(2C cu-co ) 

r f ’r f r 


(6.47) 


or as 


q r (t) - c " ?r “ r ' :|C r cos “dr t + D r sin “dr t > + + 

r 


j (os t+0 ) 
, »r> Xr>° f 1 

J (0l r S -u f 2 ) 2 +(2? r ^u r ) :i 


, r - 1,2,3, . . . ,n. (6.48) 


For convenience in further calculations, the equation above can be 
rewritten as 


-5.05 t 

{^(t)} ■ e [C r cos05 dr t + D^sinto^t] + j[$] {F q } + 

2w r 

[$] T (F n }sin(u) t+0 ) [4>] T {F n }cos(u)-t+0 ) 

4 - r r ^ Uc f r 

/(0) 2 -W S 4 +(2C 0) 05.) 2 / (iu 2 -0).^) 2 +(2? 05 05.) 2 

r f *r r f' r f r r f 


r ■ 1,2,3, . . . ,n. 


(6.49) 


where: 



ORIGINAL PAGE IS 
OF POOR QUALITY 

(f } ■ / {f* + {f T 

u 0r 7 f U 0s ' + u 0c' 


95 


(6.50) 


ij» r ■ arc tan 


6 “ arc tan 

r 


'Os 


F 0c 
ju -Jr 


2 2 
(A) -UL. 

r f 


, r * # » • jHj (6. 51) 

, r - 1,2, . . . ,n, (6.52) 


and 


0 


T m % - V r “ 1,2,3, ...,n. 


(6.53) 


To define the vectors of constants (c) and {d} , we apply the 
Initial conditions given as the vectors of initial displacements and 
initial velocities (x(0)} and (x(0)}. To apply the Initial conditions 
from the original system to a new, transformed system, we use eqs. 
(5.69) and (5.70) at t ■ 0. 


(x(0) } - («l>){q(0)} 


(6.54) 


(x(0)} - [$){q(0)>. 


From equations (6.54) and (6,55), we have 


- 1 , 


(q(0)> - W f (x(0) } 

(q(o) } - i*i" l {x(o)}. 


(6.55) 


(6.56) 

(6.57) 


In order to avoid the inversion of matrix [$] , we are using the 
orthogonality relation derived in Chapter 5 


t*rimH*] - (I! 


(6.58) 
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where [I] is the identity matrix. Postmultiplying equation above 
by we have 

t*)" 1 - (6.59) 


Substituting eq. (6.59) into eqs. (6.56) and (6.57), we have 

(q(0) } - [4>] T [m] (x(0) } (6.6Q) 

{q(0)> - [$] T [m]{x(0)>. (6.61) 

At the time t - 0, eq, (6.49) reduces to 


{4^} i {F Os }sin0_ i 


i T r Os r , 

q ( 0 ) - c + j{9 } {F } 2 2-2 2 " 

2o» Z r 0 /(w Z -w/r+(2C w U)J 2 


r f 


-* W W f, J 

r r f 


{ V {F Oc }cos0 r 


_, r - 1,2,3,.. .,n. (6.62) 


/ (w 2 -(0, 2 )^ 


(oi r ‘-w f “) +(2C r W r u) f )‘ 


Substituting eq. (6.60) into eq. (6.62), we calculate the vector of 
constants {c} 


C 


r 


M r } T ta]ji r (0) 



<y T tv 


r 


{ y T{F 0 s }slne r 

/ (w r 2 -w f 2 ) 2 +(2? r a) r 0) f ) 2 


{ y T{F 0 c }cos9 r 

/ (co r ^ f i ) 2 +(2C r a )r a3 f ) 2 


(6.63) 


Finding the first derivative of eq, (6,49), we have 
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-? u t 

Cq r (t)} - -? r w r e r r tC r cos» dr t + D r sinw dr t] + 


-? w t 

+ e r r [-C r W dr sino) dr t + D r 0) dr cosw dr t] + 


0)^{$^} T {FQ g }cos(Wjt+8 r ) 0) f {$ r } T {FQ C }sin(u)^t+8 r ) 

/7w r 2 -w f 2 ) 2 +(2? r o» r U) f ) 2 / (w r 2 -w f 2 ) 2 +(2^.w r iiJ f ) 2 


r * l,2,3,«..,n« 


(6.64) 


At time t - 0, eq. (6,64) becomes 


{q (0) } ■ •? li) C + D to , + 

H r v ' r r r r dr 


“f { V {F 0s }cos9 r 
/ (w_ 2 -<o, 2 ) 2 +(2? r a) r a) f ) : 


r f 


VV {F Qc }sln9 r 

/(a) r 2 -a) f 2 ) 2 +(2C r UJ r a) £ ) 2 


, r - 1,2,3, . . . ,n. (6.65) 


Substituting eq. (6.62) into eq. (6.65), we calculate the vector of 
constants D^, and in a similar manner, we calculate the vector of 
constants C r : 


r to 


dr 


T • MO ^ F Os }cOS0 r 

{4> r } T tm]x r (0) - 5 r W r C r - 


/ (to 2 -<o 2 ) 2 +(2? to to.) 2 


r f 


*r r f 


w,{* ) {F n >sin8 
f r 0c r 


/ (w r 2 -co f 2 ) 2 +(2? r to r to f ) 2 


, r * 1,2,3, ...,n. (6.66) 
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The vector is defined by eq. (6.63). Now we can find the comDlete 
solution vector of the transformed system (q(t)}. Finally, the dis- 
placements of masses are obtained by superposition of the response in 
each mode 

(x(t)> - m(q(t)}. (6.67) 

Very often, only the steady state response of eq, (6.49) is of 
interest to us. For this case, x r (0) - 0, and x r (0) - 0, r « 1,2,3,..., 
n, and eq. (6.49) is reduced to only the particular solution, 


q r (t) 



(♦/V + 


{°r }T{F 0s }8ln(ai f t+e r ) 

/(m i -w» 2 ) 2 +(2? 0) w y 

r f r r f 


+ 


{* r } T {F Qc }co s ( t+8^ ) 

/ (0) 2 -oj V+(2C « w f )* 

r f r r f 


1,2,3, 


9 


n. 


( 6 . 68 ) 


The vector of phase angles between the harmonic forcing-function 
vectors and the displacement vectors, 6^, r » 1,2,3, . . . ,n, is given by 
eq. (5.62). The modal amplitudes of sine and cosine forcing functions 


are given as 


«s> 


{<& r } T {F Og }sin(w f t+0 r ) 

m ■' I - n. l— 11 1— — ■ ■ ■!. — » ■ m 

r / (a) r 2 -a) f 2 ) i +(25 r to r (D f ) 2_ 


9 ^ m 1 9 ^ > 3 9 • • • , tl , 


(6.69) 


and 


1 


a 
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(Q c > 


r 


{4» r }T (F 0g } cos(w f t+6 r ) 

/(ft) 2 -w 2 ) 2 +( 2? w w ) 2 

r f r f' 


r « 1,2,3, . . . ,n. 


(6.70) 


The dynamic amplification factor in a mode is defined by the analogy 
for single-degree-of-freedom system, given in Appendix A. The vector 
of the dynamic magnification factors in modes is defined as the ratios 
of the steady-st&te modal amplitude to the vector of the static 
modal deflections (Q ) : 

* SC 

Q t 

D * , » r “ 1»2 ,3 , . . ,n. (6.71) 

9 r ( Vst 

The quantity above is usually given for a particular component of a 
harmonic forcing function. The response due to a sine forcing function 
in eq. (6.68) can be written in terms of modal amplitude: 


q r (t) " Q r sin(a) f t: + V» r " l,2,3,...,n, 


where Q r is the vector of modal amplitudes of the sine function, written 


as. 


{ V {F 0s ,sin9 c 




-, r ■ 1,2,3,. ..,n (6.72) 


(ai -u ) +(2C ft) ft)-)‘ 
r r r r i 


Observing eq. (6.34), the static displacement at the r-th mode is 


Vst ■ iT 


j(<l> r } T {F 0 } + {$ r } T {F 0g }sini|i r + {4» r } T {F 0c >cos^ r 


i 


(6.73) 
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We can expand equation (6.73) for a very general case: 


^r^st * ( 1)^2 


|-[<& r ] T {F 0 > + {$ r } T {F 0s )sintp r + {<t r } T {F 0c }cost|/ r 


r • 1,2,3,. ..,n. < 6 * 74 > 

The maximum dynamic steady-state displacement from eq. (6.68) is 


, T <#„> T tF n > 

Q r ' 77T <Q r J T <F 0 >+ «*- 

r 


/ (w. 2 -^ 2 ) 2 +(2? r 0) r W f ) 2 


r f 


r - 1,2,3, ...,n. 


(6.75) 


The dynamic magnification factor in a mode is 


CO. 


„ T St . -loo 

D r - — ^ , r * 1,2,3,.. .,n. 


(6.76) 


As a very important characteristic of the system, the damping 
in a mode can be defined by one of several methods given in Appendix A 
for a damped single-degree-of-freedom system. Amount of damping in a 
mode is usually determined by measurement for a particular structure 
or by analogy with the damping of a similar structure. To define the 
damping in a mode, eq. (6.34) can be set for a free vibration: 

q r (t) - 2C r CJ r ^ r (t) + U) r 2 q - 0, r - 1,2,3 n. (6.77) 

Letting q r (t) * A^e^, eq. (6.77) becomes 

[p 2 + 2? 0) p + (0 2 JA e pt - 
^ r r^ r r 


0. 
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Solving the equation above, we have 
p l,2 " ■ 'r r 


+ i /h„V 2 - 4w 2 , r - 1,2,3, (6,78) 


We concluded In Appendix A, that transition between the imaginary 
and real roots of the characteristic equation (eq. (6.78)) occurs when 


or 


2 , , 2 ... 2 

4C w ■ 4w , 
r r r * 


C r - i.o, 


(6.79) 


where is the ratio of the actual damping in the r-th mode to 
the critical damping in the r-th mode 

C 


’ r < C nrh * 

cr r 


(6.80) 


From eq. (6.79), the critical damping in the r-th mode is 


(C ) ■ 2w . 

cr r r 


(6.81) 


The actual damping in the r-th mode is 


C - 2? a) , 
r *r r’ 


or, for a very general case, 


(C cr ) r m 2ta r’ ** * 1,2 ,3, ... ,n, 


and 


C - I t, oj , r 
r r r 


1,2,3, .. « , n . 


(6.82) 


(6.83) 


(6.84) 
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Again, by the analogy given in Appendix A, the vibration in a mode 
is conventionally classified with respect to the degree of damping 
as 

1, Critically damped when * 1.0, 

2, Overdamped when ? r > 1.0, 

3, Underdamped when £ r < 1.0. 

In previous derivations we assumed that we had an underdamped 
case or oscillatory motion, which is the only case of interest in 
engineering application. 

Once more, the same question arises regarding how many modes 

we must take in order to have the solution accurate enough for an 

engineering analysis# Since, for the dynamic cases, inertial 

forces of the system are very important, we usually want to know 

* 

how well we model our system dynamically,; Since we usually discretize 
the system in modal analysis, a very good indication of amount of 
mass taken into consideration is "the effective modal mass." For 
the r-th mode, this quantity is defined as 

(1} T m {*}{* } T [mJ {1} 

CK ef f) r - — I — ■ - ■> r = 1,2,3, ...,n, (6.85) 


or, in a summation form. 


(M cc ) = 
eff r 


n 


£ m. <J> 
i»i 1 ir 


n 

E m $ ' 

. , i ir 
i=l 


-, r = 1,2,3, ...,n. (6.86) 


X 

Where m^. is the mass at the i-th point, {l } is a row matrix of ones, 

and M is the modal mass in the r-th mode defined by 
r J 
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2 

M * 2 m , $ . , r “ 1|2 (6.87) 

r ±ml i ir 

A good Indication that the inertia forces have been accurately 
measured in modal analysis is usually that the effective mass is 
bt least 95% of the actual mass of the system. Since the modal 
parameters for higher modes cannot be calculated correctly, inclusion 
of higher modes usually does little to. increase the accuracy of 
the calculation. For structural systems with homogenous material 
characteristics, the first several modes are most often enough to 
calculate the dynamic response of the system. 
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Chapter 7 

MODAL ANALYSIS PERFORMED ON FIVE BLADE MODELS USING A 
NICOLET 6 60 A DUAL-CHANNEL FFT ANALYZER 

7,1. Introduction 

In the previous chapter, we discussed the concepts which have 
to be considered during modal testing using FFT analyzers. In this 
chapter, we are going to give a brief description of how to conduct 
the modal testing properly using a Nicolet 660A Dual-Channel FFT 
Analyzer when Impulsive excitation is chosen. Selection and 
calibration of the impact hammer will also be discussed, and the 
results of testing performed on five models will be presented, 

7.2, An Experimental Technique 

Impulsive testing utilizes short-duration transient force Inputs 
(with corresponding broadband spectra) to excite all frequencies in 
the structure simultaneously for a chosen frequency range. Modern 
measurement Instrumentation (such as the Nicolet 660A) is fully 
capable of performing the modal analysis when the impulsive 
excitation is used. 

Impulsive testing is normally conducted using some form of an 
instrumented hammer. A force transducer is mounted on the head of 
the hammer, and measures the force input to the structure. An 
accelerometer is used to measure the response of the structure. By 
far the most popular technique for impulsive excitation is to mount 
an accelerometer at a fixed location and excite the structure at a 
multiplicity of locations, using a hammer with a force transducer on its 
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face co Impulse Che syscem at many locations. The daCa are analyzed 
wich a dual-channel FFT analyzer, using Che Impulse force as a cransienC 
capcure crigger condicion. These two transient time -histories are 
Fourier-transformed to yield the Input and output spectra. The 
resultant ratio of the output and Input spectra is the desired frequency- 
response function. 

Usually, the results of several transient excitations are 

ensemble-averaged. Then the signals are Fourier-transformed and the 

transfer function is displayed on the screen of the analyzer. The 

Impulse hammer (excitation system) with sensitivity S (mV/lb) is 

a 

attached to channel A of the FFT analyzer. The accelerometer (response 

transducer), with a sensitivity of S OaV/g), is attached to channel 3 

b 

of the analyzer. Based upon a preliminary examination of the structure, 
the input attenuator of channel A is set to a full-scale value of F , 
and the input attenuator of channel B is set to a level of F^. 

The levels of F^ and F^ can only be obtained by trial measurements 
on the structure. 

The structure to be examined is marked &t M spaclal locations. 

The response transducer (.accelerometer) is fastened to the structure 
at one location. This transducer should be located in a local "hot 
spot" so that it will be capable of responding to all inodes. Starting 
from the response-transducer location, the structure is impacted a 
number of times at each spatial location. Each spatial location is 
impacted until the coherence function (a function which indicates 
the validity of the transfer function) approaches unity across the 
selected frequency-analysis range. The first impact will always produce 
a coherence function equal to unity. It is essential that each 


location is Impacted more than once, so that the coherence function 
may be used as a valid indicator of the quality of acquired data. 
The first impact point is of particular importance. 
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A complete description of a mode contains four pieces of infor- 
mation. At the n-th mode: 

1. Natural frequency of the n-th mode (f R ) , 

2. Amplification factor of the n-th mode (Q^), 

3. "Weight" of the n-th mode; an indication of how much 

structural mass is in motion in the n-th mode (W 5, 

n 

4. Mqtion of the m-th. location in the n-th. mode 

n >. 

m,n 

Within a given frequency range, N resonances may exist. Hence, 
n ranges from 1 to N. If M locations within the structure are 
measured, then m ranges from 1 to M. Each of the N mode shapes may 
be normalized so that the maximum displacement (.for example, at the 
p-th location) is equal to one unit. The choice of mode shape 
normalization affects the moda,„ weight W fl . The modal weight is 


defined by the summation 


M 

w - t w (<t> r 

n m-l n m > n 


(7,1) 


where W is the physical weight at the m-th location, 
m 

At the first impact point, both input and response readings 
are made at the same spacial location. The input/output relation- 
ships associated with this spacial location are used as a guide to 
interpret all subsequent measurements. The resonant frequencies 
of a structure are independent of the spacial location from which 
they are measured. The resonant frequencies of N modes are detected 
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by identifying the N peaks in the imaginary part of the displayed 
transfer function. The quality factor Q's of the N modes are 
measured using information from the displayed real part of transfer 
function, and are defined as 


where 



a bn 


<f /f.) „ 

a bn 


+ 1 
- 1 


(7.2) 


f - frequency above resonance, where the real part 
a 

of the transfer function reaches a peak, 

and 

f^ ■ frequency below resonance, where the real part of the 
transfer function reaches a peak of opposite sign. 

We showed in section 7,4.3 that the quality factor in a mode -is 
directly related to the amount of damping in a mode. The first 
measured transfer function, identifies the parameter f and 0 over the. 
range of 1 < n < N. 

Mode shapes can be estimated directly from frequency-response func- 
tions. The mode shape is estimated by measuring quadrature (imaginary) 
response values from all of the measured transfer functions at each 
resonant frequency. At each of these frequencies, the quadrature 
components of the response functions are distributed in proportion 
to the displacement at each of the measured locations throughout the 
structure. That is, they are distributed in an approximation of the 
mode shape associated with that resonance. 

The imaginary values from the driving-point transfer function at 
each of the N resonances should be measured and retained to normalize 
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subsequent data. The N quadra:ure points, XMAC^^Cf^) , must be 

retained. The driving-point location, a, is subsequently moved 

through the sequence of predefined spatial locations. At each of the 

M locations, the Imaginary part of the transfer function is measured 

at the N preselected resonant frequencies, The quadrature values, 

XMAG , (f ) , are measured and retained for all values between 1 and 
a,b a * 

M and all values of n between 1 and N. 

If all measurements can be acquired without changing the input 
attenuators of channel A or channel B, then a single scale factor is 
sufficient to relate physical and electrical units. The constant, 
k c , may be retained where 


k 


c 


- S b F . /S aV 


(7.3) 


The coefficient, k c , relates the voltage gain measured by the 
analyzer to physical measurements in the structure. If the imaginary 
part of* the transfer function, IMAG fl j j (f n ), is read at the n-th 
resonant frequency, it is related to the physical parameters of the 
system by the equation 

k 0 <|> <J>, 

IMAG (f ) - ■ - - (7.4) 

a ,b n w 

' n 

At any frequency, £ r , the stored quadrature components, IMAG a ^(f^), 
will have a local minimum occurring at some special location a ■ p. 

The components IMAG^ ^(f^) are of particular importance. These are 
used to compute the "weight" of each mode and to normalize the 
retained quadrature components from each mode to arrive at the spatial 
vectors <f> : 
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k A IMM bb«n> 

M *»■ ■ . ■■■ ■■ bun % 

IMAG* .(f) 

P|b n 


♦ 

a,n 


^-■b^ 
““p.b'V ' 


(7.5) 


(7.6) 


These shape functions are approximations, because the measured shape 
(in response to a single-point excitation) is actually ?, linear summation 
of all of the mode shapes of the system. If the resonances are broadly 
spaced In frequency and high M Q" in nature (Q « 5 or greater) , then 
the interaction between the modes is minimal and the distortion of any 
given mode shape by contribution from other modes is not significant. 

7.3, Calibration of an Impact Hammer 
An impact hammer of 12.5 [lb] equipped with a light Teflon tip 
and the Kistler Instrument Quartz Force Transducer Model 901A-903 
were selected, and the force transducer output was fed directly to 
channel A of a Nicolet 660A Dual -Channel FFT analyzer. In Fig. 7.1, 
the FFT analyzer and impact hammer kit (impact tips and weight 
extenders) are shown. Preliminary testing on models showed that an 
impact hammer with a light Teflon tip and without weight extenders 
was enough to excite all frequencies in a frequency range of 0-2 KHz 
for all five models. A freely supported mass of 10.437 [lb] (steel 
cylinder) was impacted. An Endevco Piezoelectronic Accelerometer, 

Model 2233, was mounted at the center of one face of a steel cylinder 
of an impact hammer and the Teflon tip at the center of the opposite 
side was used as an impact surface. (The mass was impacted through 
.the center of gravity.) A piezoelectric accelerometer is a self- 
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generating transducer which generates an electrical output signal 
that is proportional to its acceleration. An Endevco Input Amplifier, 
Model 2614B, with an amplification of 10 was used to amplify the 
signal generated by the accelerometer. This signal was fed to channel 
B of the FFT analyzer. The whole calibration was performed in a 
frequency domain of 0-2 KHz. 

An ensemble average of a number of impacts was done, and the 
calibration factor was calculated in the following manner: 

The unitless constant k was calculated as: 


where 


c " A (to) 


. 6950.65 

c 486.73 


14.296 


k - 14.296 
c 


(7.6) 


F(to) is the average value of the Fourier -transformed force 
signal in a frequency range of 0-2 KHz, 

A(to) is the average value of the Fourier-transformed acceler- 
ation signal in a frequency range of 0-2 KHz. 

The mass of an impacted cylinder was calculated as: 


m *— ■ 0.0270 [lbsec^/in], 
c q 

The sensitivity of the accelerometer was 


(7.7) 


E a « 37.5 ImV/g]. 


(7.8) 


The total sensitivity of an accelerometer in engineering units was 
found by 


{ 
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S - 


E /2/2 

a ' 


a Ampl. Factor 


[V/in/sec ]. 

/r 


(7.9) 


The expression above was multiplied by y, since the signals were 
displayed on the screen of FFT by their RMS values. 

0 37.5x10" 3 x /T /2 


S - 0.686x10 
a 


10 


-5 


in/sec 


(7.10) 


This was the total sensitivity of the accelerometer and in-line 
amplifier with a voltage gain of 10. 

The sensitivity of the force transducer was not known. The 
force transducer system was calculated using eq. (b) derived in 
Appendix B. The two total sensitivities are related as follows: 


F(u>) , _ m c S f 

A(w) c S 

a 


(7.11) 


where 


and 


m c is the mass of an impacted cylinder 


is the total sensitivity of the force-transducer impact- 


hammer system. 

Thus, from eq. (7.11), we obtain 

k S g 


c a 


[mV/lb ] 


S f " 


w f G 

14.296 x 37.5 x 10~ 3 x *ff/2 
10.437x10 

/r 


(7.12) 


Again, eq. (7.12) is multiplied by y to get the values on the screen 


i 
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in engineering units (zero to peak). 

S f - 0.363 x lCf 2 [V/lb] (7.13) 

This technique has the added advantage of being able to calibrate 
an accelerometer and impact hammer as a matched set. If the 
accelerometer used to measure the response of the suspended mass is 
the same accelerometer that is used during the modal test, then any 
error in the accelerometer measurement for the suspended mass will 
automatically be introduced into the impact hammer calibration. 

This error will then be eliminated when the ratio of acceleration to 
force is calculated, since it will be in both the numerator and 
denominator. Therefore, the accuracy of the calibration of the 
matched set of transducers depends on the accuracy to which the 
calibration mass is known. A very good low-frequency calibration 
can be obtained using this technique. The high-frequency limit is 
determined by the frequency response of the accelerometer and the 
mass. Figures 7.2 and 7,3 show the setup for the Impact hammer 
cal ibration. 

7 .4 . Geometrical and Material Properties of Blade Models 

The modal testing was performed on five blade models. Three 
aluminum blade models had NACA four -digit symmetrical airfoils as 
their cross-sectional ar*.«,:». The other two blade models were of 
aluminum and steel, w;td had identical rectangular cross-sectional 
areas and comparable lengths. The problem was to relate the damping 

of the geometrically identical models made of different materials, and 

# 

to observe the damping of the models made of the same material 



Figure 7.2. Schematic of test setup for calibration of the Impact hammer. 
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figure 7.3. Test setup for calibration of the impact hammer 
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(aluminum), but having different cross-sectional areas. 

As a very important concept, damping for the same family of 
models could be much more easily handled if we could get it in a 
scaled form with respect to the chord and the thickness. These 
blade models are shown in Fig, 7.4. Geometrical and material 
properties of these models are given in Table 7.1. The sixth model 
was an aluminum curved blade, as shown in Fig. 7.5. The objective of 
this test was just to measure natural frequencies. A different 
setup was used for testing of the sixth model. 

7.5. The Modal-Testing Procedure 

The modal testing was performed on the five blade models 
described earlier by the Nicolet 660A Dual-Channel FFT Analyzer using 
impulsive excitation. The accelerometer was bonded beneath the tip 
of the blade (using adhesive). The "target” points were bonded to 
the upper surface of the blade, along the elastic line (a locus of 
shear centers) . Only four points equally spaced counting from the 
tip were used as impacting points. Since the data were analyzed 
over a 2 KHz range, all frequencies should be excited in that range. 
To accomplish this, extreme care was taken to locate the target 
points properly. The nodal points for the first five modes were 
defined and good arrangement of the impact points was done to avoid 
their location at any of the nodal points. Figure 7.6 shows the 
arrangement of the nodal points of the blades for the first five 
modes. For this type of model, the first 5-6 natural frequencies 
are within the chosen frequency range of 2 KHz. 
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Geometrical and material properties of blade models used for modal testing. 
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The accelerometer’s signal was conditioned by an Endevco Input 
Amplifier, Model 2616B, with a voltage gain of 10. During the 
calibration procedure, we defined the total sensitivity as 


S « 0.686 x 10 
a 


-5 


in/sec 

Impacting the targets wit 


JL 

lb 


The blade models were excited by 
i an impact hammer, shown in Fig. 7.1(b). 

A force transducer, Kistler Instrument Quartz, Model 901A-903, 
was mounted to the machined face of the hammer with a threaded stud. 

A plastic tip (Teflon) was mounted to the force transducer to serve 
as the striking surface. The reason for Its mounting was discussed 
in the previous chapter. During the calibration procedure, we defined 
the total sensitivity of the force transducer and Impact: hammer 
system as ■ 0.363 x 10 " 

We discussed earlier the importance of ensemble-averaging both 
signals in the time domain. When the signals recovered from noise, 
they were Fourier -transformed and the transfer function was calculated. 
To obtain above both transient signals in the time domain, they 
were captured and stored. The best way to do this is by automatic 
capturing of the signals. Working with this analyzer, we set the 
trigger on channel A (force signal) to activate the system (capturing 
and storing) of both signals at the of the maximum amplitude. . 

Also, the type of averaging was chosen so that besides the automatic 
triggering, capturing and storing the averaging was done automatically. 
The blade was struck a number of times until the satisfactory 
coherence function was observed on the screen of the FFT analyzer. 

The stable data were considered when the coherence function exceeded 
a value of 0.9 around each resonance. This was an indication of 
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valid data free from nonlinear effects and noise contamination. 

Since the trigger was very sensitive on the level of the force signal, 
enough time was allowed between two impacts for the transient 
vibration of the model to die out completely. By our observation, 
this period was usually to two minutes. The whole process of 
averaging and the capturing of signals was observed on the screen. 
When good coherence was reached, the transfer function was calculated 
and displayed on the screen. Also, this function could be described 
in all three modes: Bode plot, Nyquist plot, and real and imaginary 

part vs. frequency. Then, the picture from the screen was reproduced 
on the paper by the plotter attached to the FFT analyzer. The same 
procedure was repeated for all target points of the blade. The setup 
for modal testing used is shown in Fig. 7.7. Figure 7.8 shows all 
components of the setup ready for testing. This test was performed 
In the structural area of the Department of Aerospace Engineering 
of The Pennsylvania State University during February and March, 1982. 

7.6. The Results Obtained by Modal Testing 

The objective of this test was to analyze the natural frequencies 
and modal damping of the blade models. Ths mode shapes were not 
analyzed, since a very good conditioning of both signals was required. 
The accuracy in measurement of mode shapes directly depends on the 
accuracy of the determined relationship between electrical and 
engineering units. For accurate measurement of frequencies and 
damping, we do not need such an accurate conditioning of the signals. 

Frequencies and damping were calculated from the point transfer- 
function of the blade. The point transfer-function is calculated 




Figure 7.7. Schematic of the setup for modal testing using the Nlcolet 660A Dual Channel 
Analyzer. 
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The setup for modal testing; (a) the blade model of the 
test stand, (b) Che components of the setup. 
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for the station where the input and output signals are taken from 
the same point on the blade. In our ease, this point was the first 
target point at the tip (the target and accelerometer were aligned 
at the same axis). Since the natural frequencies and damping are 
the global property of the linear structure, any arrangement of the 
impact point and accelerometer should give the same natural frequencies 
and damping. The measurements and plots were made for four different 
locations of impact points, and good agreement was observed. 

When good coherence was reached, the transfer function was dis- 
played on the screen of the FFT analyzer, and the frequencies were 
read accurately, hut using the cursor line, as digital numbers. 

Figure 7.2 shows the moments of impacting with the blade. Figure 
7.9(a) shows impacting with target No. 1, and Fig. 7.9(b) shows 
impacting with target No. 4. Figure 7.10 presents the test setup 
for modal testing of a steel blade. Since we do not have enough 
space to show all graphs reproduced from the analyzer screen, only 
those graphs typical of all blade models will be shown 
including an aluminum blade model having NACA 0015 as its cross 
section, and a steel blade model having a rectangular cross section. 

The complete results for all five models are given in Table 
7.2. They will be used in the following chapter in the numerical 
model of the dynamic response of the blades. At this point, we 
cannot make any conclusion concerning the data measured. A preliminary 
check shows that at least we have good agreement regarding natural 
frequencies compared with the natural frequencies of continuous 
beam models. Also, we observe that the damping is concentrated 
in the first three or four modes. 
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Figure 7.11. The magnitude and phase of the point transfer function of a model having an NACA 0021 
cross section (peaks of the magnitude correspond to the model's natural frequencies). 
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Figure 7.13. The coherence function related to Input/output used for determining the point transfer 
function of a model having an NACA 0015 airfoil cross section. 
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Figure 7.14 . The magnitude and phase of the point transfer function of a model having an NACA 0015 
cross section (peaks of the magnitude correspond to the model's natural frequencies). 
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Figure 7.15. The real and Imaginary parts of the point transfer function of a model having an NACA 
0015 airfoil cross section*. 




Figure 7.16. The coherence function related to the Input/output used for determining the point 
transfer function of a model having an NACA 0015 airfoil cross section. 
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Fieure 7 19. The coherence function related to the Input/output used for determining the point 
transfer function of a model having an NACA 0009 airfoil cross section. 
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Figure 7.20. The magnitude and phase of the point transfer function of an aluminum model having 
rectangular cross section (peaks of the magnitude correspond to the model's natura 1 
frequencies) . 
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Figure 7.23. The real and Imaginary parts of the point transfer function of a steel model having 
a rectangular cross section. 




Table 7.2. Natural frequencies and nodal danplng ratios obtained by nodal testing. 
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(a) A curved blade model at the test stand, (b) the 
test setup for measuring frequencies of a curved blade 


Figure 7.24. 
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Chapter 8 

NUMERICAL MODELINC AND RESULTS OF THE DYNAMIC FLEXURAL 
RESPONSE OF UNIFORM PROPELLER BLADES 

8.1. Introduction 

This chapter 1* a summary and discussion of the results and 
findings of this study. A small computer code was developed to 
perform a complete modal analysis of all five blade models, including 
the dynamic response of the models as undamped and damped multi- 
degree-of-f reedom systems. The data obtained experimentally by 
using the FFT analyzer, which were analyzed in Chapter 7, will 
used for a damped case. Some important results of numerical modal 
analysis will be shown and discussed. A very brief description of 
the computer program will be given in the last section. A flowchart 
of the computer program is given in Appendix C. 

8.2. Formulation of the Equation of Motion 

Since we are going to model the blades as undamped multl-degree- 
of-freedom f ree-vibratlon cases, the formulations of the (ml and [k] 
matrices are necessary. In section 5.4, we gave a detailed description 
of methods for deriving the mass matrix [m]. We found it convenient 
to assemble this matrix as a diagonal matrix having only diagonal 
terms and lumped masses. A schematic representation of this model 
is given in Fig. 8.1. 

The formulation of the stiffness matrix [k] is done in the 
following manner. First, the matrix of flexibility influence 
coefficients is found; then, by inversion of this matrix, [k] is 
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Figure 8.1. (a) Schematic of the blade model, (b) masses lumped at 

each ncde of the model segment, (c) final arrangement of 
the lumped masses of the blade model. 
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computed. Both of them have to be symmetric matrices. The matrix (k) 
can also be calculated by using the direct stiffness method. The 
computer program Is capable of handling the equation system when 
we have a full-mass matrix (m) . The program la run for two cases: 
when the order of these matrices Is 8x8, and when It la 10x10. The 
program la capable of handling larger arrays. 

8.3. Eigenvalues, Natural Frequencies, and Eigenvectors 
Since the accurate calculation of these parameters la very 
important for the dynamic analysis of the models, extreme care and 
rigorousness was taken In the process of setting up the algorithm 
of computer program. The governing equations of motion can be written 
as 

[m] (x(t) } ♦ (k]{x(t)} ■ {0}. (8.1) 

The Jacobi method Is used to calculate all eigenvalues and corresponding 
eigenvectors of the generalized elgenproblem 

(kj M>] - w 2 [m] [$]. (8.2) 

The program basically transforms the Initial system of equation (8.1) 
into a new system having only the diagonal elements In both matrices. 

This Is a so-called modal system of equations, using the modal 
mass matrix [MJ and modal stiffness matrix [KJ Instead of matrices 
[m] and [k]. In Chapter S, we showed that the following relationship 
exists : 




[0J tmj [0] - [Hj, 


(8.3) 
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and 

»J CkJ W - CKJ (8.4) 

where the matrix (<f>] is called the modal matrix or array of eigen- 
vectors of the system. Also, in Chapter 5 we showed that we could 
scale the eigenvector matrix In a certain manner to get the Identity 
matrix [I] instead of the matrix [MJ . We established the following 
relations : 


[♦] [m] [$] - CD, 

(8.5) 

[$j m w - cA, 

(8.6) 


where ($] is the matrix of scaled eigenvectors. We showed also 
that It Is much more convenient to work with unit masses, where we 
could avoid some computational expenses during the dynamic analysis 
of the system, especially for imposing the initial conditions of 
the system. We will be working with the scaled system of equations 
throughout the rest of the program. In order to visualize the 
deformation pattern, it is convenient to scale eigenvectors in a 
similar manner as to how unit vectors are defined. Then they have 
values only between zero and one, and it is easy to see the relative 
motion of certain points of the system in different modes. The 
natural frequencies are also calculated when the blades are modelled 
as continuous beams (Bernoulli-Euler Model). In Table 8.1, the 
frequencies are given for all five models for the first eight modes. 
Comparison is usually satisfactory among the three approaches used 
up to the fifth mode. A continuous system and lumped parameter 
system are used as numerical approaches to identify the first eight 


Table 8.1. Predicted and Measured natural frequencies for five blade Models. 
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frequencies. Only the first five frequencies are measured by the 
FFT analyzer (since the frequency range chosen was 0-2 KHz). 

In Fig. 8.2, the deformation pattern is given for the first five 
modes (greatly exaggerated) . Also, pictorial representation of 
every element of the modal matrix is given . This is given only for 
five modes, although the pattern is the same for any number of modes. 
In Table 8.2, the normalized mode shapes are given for blades having 
NACA 0015 airfoil cross sections. The system of notation corresponds 
exactly to the pattern given in Fig. 8.2. 

8.4. Identification of Modal Parameters Including Damping 

In a previous section, we showed the natural frequencies and 
modal shapes for some of the blade models we studied. In this 
section, we are going to introduce the damping measured and analyzed 
in Chapter 7, then perform a complete modal analysis of blade models 
modeled as damped multi-degree-of-freedom systems. We assume that 
the natural frequencies and modal shapes are not seriously affected 
by introducing the damping, which is valid for systems having light 
damping. As we described earlier, the damping ratios are measured 
in the first five modes. Measured data cannot be used directly in 
an analytical model without modification. 

In Chapter 6, it was shown that we have to impose some 
conditions under which the modal matrix decouples the equation of 
motion for the damped multi-degree-of-freedom system. The damping 
forces could be decoupled only if the damping matrix is a linear 
combination of the mass matrix [m] and stiffness matrix [k] . The 
damping ratios are measured only in the first five modes. In the 









Mode 
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analytical model, the damping ratios in the rest of the modes are 
taken to be zero, and the least-squares approximation straight line 
fits through these points, starting at the origin. The damping 
ratios modeled in this manner satisfy eq, (6,13), which is rewritten 
here as 

{♦jJlcJWj} - (8.7) 

where ^ is the modal damping ratio and 5 is the Kronecker delta 

“ 1* for i ■ j, ” 0, for i / j). Using the transformation 

matrix of scaled eigenvectors [$] and eqs. (8.5) and (8.6), the 
governing equation of free vibration of an undamped multi-degree-o/f- 
freedom system becomes 

[l]{q(t)}+ t2WCJ{q(t)J + rU)^j{q(t)> - {0>, . ((8.8) 

where the relation between original and modal displacements is given 
by 

(x(t)} - [<J>] {q(t) }. (8.9) 

We have defined the necessary parameters for uncoupling the 
equations of motion. Since we have imposed some conditions on the 
damping ratios obtained by measurement in order to decouple the 
damping forces, the related parameters such as the amplification 
factor, critical damping, damping in a mode, and damped natural 
frequencies have new values. These parameters are given for all 
five models in the first 8 modes in Table 8.3. Since we do not 
have much interest in a free vibration, most of these parameters 
will be used in the following section, where we are going to 



Table 8.3. The nodal parameters Identified for the blades modeled as damped MDOF syst 
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consider Che dynamic response of the blade models modeled as damped 
and undamped multi-degree-of -freedom systems, 

8.5. The Dynamic Response of Blade Models Modeled by Modal Analysis 
as Undamped and Damped Multl-Degree-of-Freedom Systems 
In Chapter 6, we defined analytically the dynamic response of 
blade models modeled as multi-degree-of-freedom systems. The 
excitation force .was assumed to be one harmonic of the Fourier 
series. The equation of motion for an n-degree-of -freedom system 
was given as 

[m] {x(t) } + [c]{x(t)} + [k] (x(t) } - (F(t) > - y {F q } + 
j(w ,t+\p ) 

+ {F 0r >e 1 , 1 « 1,2,3 n. (8.10) 

Quantities {F^} and are defined by 

{F 0r } ' ^ F o. }J + < F 0c }4 

and 

■ arctan 

The equations above are written assuming that the vector of the 
excitation forces can be written as 

(F(t)} - y{F Q } + {F 0g }sinw f t + {F 0c }cosU) f t. (8.13) 

T 

Premultiplying eq. (8.10) by I4>J and using linear transformation 
(x(t)} * W]{q(t)-}, eq. (8.10) becomes 


<F 0c>l 


'Vi 


y i m 1,2 ,3, ... ,n. 


( 8 . 12 ) 
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[I]{q(t)> + r2C0U(i(t)} + rw 2 -J{q(t)} - |{>t] T {F 0 } + 

I 

+ [<t>] {F 0r )e r , i - 1,2,3 n. (8.14) 

This is a so-called modal equation which decouples the system of 
equations into n Independent equations. Then we showed how to obtain 
the homogenous and particular solution of eq. (8.14). Finally, the 
displacement of the system was transformed into the original coordinate 
system by transformation (x(t)} ■ [<&]{q(t)}. 

The computer code has been developed for a simpler mode 
(regarding the excitation forces) than the analytical model developed 
in Chapter 6. For the test program, it is assumed that the 
excitation force is of the type F - sinto^t, and that it acts at 
the tip of the blade. Also, the dynamic response (steady-state) is 
defined when the blades are modeled as a continuous system. The 
main idea in this study was to compare the response obtained by 
two approaches. Also, the initial conditions are assumed to be zero. 
The procedure of imposing initial conditions was given in Chapter 6. 
This could be done easily by carrying out the calculation in a matrix 
notation. The procedure for calculating the dynamic response is 
implemented in a computer cade for cases when the blade models are 
modeled either as undamped or damped multi-degree-of-freedom systems. 
For an undamped case with a given excitation force and initial 
conditions equal zero, eq. (8.14) is reduced to 

[I] (q(t) } + tw 2 j(q(t)}» [<t] T {F 0g }sin(W f t) , 


(8,15) 
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where the term on the right-hand side represents the scaled modal 
force or the excitation force in a new transformed system. The 
steady-state response of the system given by the equation above for 
the r-th mode is given as 

W r }I{F 0„ } 

q (t) 5 ?p-sin(U>t), r ■ 1,2,3, ... ,n. (8.16) 

r 0) r 
r f 

The final displacements of masses are obtained by superposition of 
the response in each mode: 

(x(t)} - [$Hq(t)}. (8.17) 

For the damped case with a given excitation force and initial 
conditions equal zero, eq. (8.14) reduced to 
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Equation (8.20) can be written as 

, t . W r> T{ *o.>* to <V-V 

/(W 2 -W 2 ) 2 +(2 1 w w ) 2 
r f r r i 


( 8 . 22 ) 


The same procedure can be applied for any of n decoupled equations 
of motion. The final displacement of lumped masses is as: 


{x(t) } * (<!>] (q(t) } . 


(8.23) 


8.6. Predicted Results for the Dynamic Response of Blades 
In section 8.5, we gave the outline of the model used in a 
computer program. We decoupled the equation of motion, getting the 
scaled modal equation of motion. The coefficients and excitation 
forces of eq. (8.18) are given for all eight modes in Table 8.4. As 
we have already stated, the frequency of the excitation force lies 
between the first and second natural frequencies of the blade models. 

The steady-state response is predicted by two approaches: 
when the system is modeled as a continuous or as a lumped parameter 
system. The response of the tip point along the first half of the 
time period of the excitation force (at 10 instants of time) is 
given in Table 8.5. A good agreement is obtained by the two 
different approaches over the entire interval of time considered. 

To obtain an insight on the participation of damping in a 
dynamic response and relate it to the dynamic response of an undamped 
system, responses in both approaches are calculated along the two 
time periods of the excitation force. In Table 8,6, the responses 
are given by the two approaches when the excitation force reaches 



Table 8. A. The coefficients of the decoupled equation of motion of the blade models. 


original page is 

OF POOR QUALfTY 


1S8 


4 

i 

i j 

< 

" ? 
sec 

M *3 " 
S ! 

m P* 4 

i o c 

! 

• *5 
! 

ii 

^1 

V 

PS 

o* 

C 

"3 

s 

PB 

• 

O 

*3 

X 

mr 

o* 

FT 

3 

mA 

IN 

o 

% 

5 

IN 

o 

- 

SC 

0 

1 
c 

© % 
*3 *3 

^ s 

m in 

o’ o’ 

1 n 

U arj 

sc 1 

- mm 

1 ! i 
• • 

D o 

K N 

5 

3 O 

To 

I 

C 

0 

1 

o’ 

"*0 U 
i 

o’ 

ia 

" ri 

r U 

4 Cl 

To 

*3 

IN 

em 

IP 

o' 

j 

< 

pa + 
D O C 

1 M 1 

21 
5 O* C 

»% 

? 
4 IN 

> o* 

o 

M 

IN 

mm 

«4 

o 

0 

1 
e 

mm 

o 

mm 

M 

m*. 

e 

o 

A 

mrn 

d 

'© 

*3 

IN 

PA 

o* 

P 

b 

M 

HA 

« 

IN 

o’ 

e 

i 

o* 

'S % 

ii 
• • 

o o 

1 1 
* mm 

o' o’ 

o 1 

m mm i 

Li i 

d* o 

rt 

Vi 

D O* 

% 

P 

5 

IN 

o’ 

0 

1 
mrn 

d 

'Tr 

mm 

| 

a ; 

PA 1 

• 

o 

2 

P IN 

P «N 
9 O 

p 

b 

*3 

5 

o 

s | 
< 

% m c 

m *P m 

X • 

a s 
• • 
DOC 

~o 

- -m 

< M 

?« 

Li£ 

i 

0 

mm 

1 
PA 
mm 

o 

0 

1 

mrn 

o 

•a 

T5 

IN 

Pm 

o’ 

% 

mrn 

* 

43 

mm 

o 

IN 

o 

mm 

s 

•4 

o’ 

l 

o’ 

0 

mm 

1 
o* 

■S o 

mm mm 

| 3 

S rn 
• • 

o o 

•at 

mm mrn 

a 3 

o’ o’ 

>a " 

m mm 

l* \ 

D* O 

b 

1 
N HA 

3* O 

\ 

p 

2 

Pm 

o' 

o 

-m 

2 

d 

~o~< 

*m 

£ 2 

PA 1 

o’ 

s 3 
5 a 

3 o’ 

•o 

3 

p 

o’ 

s i 

c 

* 2 
DOC 
s • 

5 i 

5 r» C 

m mm 

i o c 

a 
» 1 
r 1 
l-l 

w 

LI 

«p 

■o 

p4 

M 

s 

mm 

o* 

o 

j 

o’ 

Is 

*3 

* 

♦N 

o’ 

o 

• 

«p 

o 

~c 

•m 

H 

5 

i 

a 

c 

i 

d 

% % 

N N 

3 a 

mm mrn 

o’ o’ 

1 

a 

ml 

Cl 

a 

B 

rW 

yu 

1 

0 

1 

o’ 

- ( 

IN 

mrn 

6 

9 'o 

sl 
s ' ?■ 


1 j 
< 

m mm m 

■ M 

« : 
m mm m 

D O C 

a 

P 4 

rl 

rl 

u 

a 

0 

1 

mrn 

o 

•o 

M 

HA 

4> 

IN 

o 

i 

mrn 

M 

mm 

PA 

o 

•* 

mrn 

8 

O 

V 

3 

s 

o’ 

0 

1 

d 

o -o 

mrn rnm 

1 2 
rnm mm 

o’ o 

rnm mm 

X X 

s a 

o’ o’ 

F*A « 

DO* 

m mrn 

§£ i 

- mrn 

D* O* 

— mm 

5 5 

- -P 

D O* 

To 

1 

PA 

o’ 

0 

1 
mrn 

O* 

mrn 

J , 

IN 

o’ 

3*0 

i 

P IN 

s' o* 

p 

p* 

3 

o’ 

S j 

c 

ll 

m mrn m 
D ©* C 

5 O 

ji 

x 

> o 

ha 

•o 

•m 

4 

IN 

PA 

o* 

0 

1 

o* 

t. 

i 

o’ 

'll 

5 

o 

0 

mm 

N 

1 

o’ 

3T" 

b> 

mrn 

M 

5 

o’ 

0 

1 

mm 

o’ 

% ~o 

PC 

H H 

S 3 

IN mrn 

o’ o’ 

%% 
mrn mm 

2 5 

o’ o’ 

m mm 

Is : 

m 

D O* 

A 

3 O 
•m mrn 

: r 

# "p 

m N 

D O* 

HA 

*o 

P 

<-> 

s 

o’ 

0 

1 
o’ 

3 

N 

rnm 

o’ 

9 O 

■m mm 

5l 

P -H 

• • 

3 O 

HA 

to 

1 

Pm 

o’ 

1 1 
p 

c 

p. 

DOC 

i IS 

m p 

D o* c 

> 

m mm 

« M 

5 

A « 

Li 

*© 

H 

X 

pa 

o 

0 

1 
p^ 

o 

O 

P 

mm 

pa 

o* 

Pm 

o 

mrn 

M 

«N 

HA 

o* 

a 

PA 

o 

*3 

*p 

IN 

o* 

0 

1 

P4 

o’ 

o I 

mm rnm 

5 3 

HA ^ 

o* o' 

*& • 
II 

o o 

D O ^ 

s? i 

m N « 

• • 

D O 

FT 

■m mm 

S 4 
3 

* *f 

1 

o 

pB 

| 

o’ 



1 

F 
* • 

i 

• j 

< 

> 2 
i zl 

m mm m 

d’ o’ c 

•3. 

m mm 

« P 

2 

- <* 

• « 

D O 

IN 

• O 

P 

00 

IN 

mrn 

o* 

0 

rnm 

1 
o* 

p* 

-9 

fN 

o 

HA 

o 

s 

«p 

o* 

M 

mrn 

R 

o’ 

m 

o 

3 

3 

o’ 

0 

1 

mrn 

o’ 

HA 

o 

mm 

s « 

o’ o 

% o 

mrn rnm 

H M 

5 s 

o* o’ 

H 

D 

m 

s ; 

55 ; 

i 2 

D p. 
ml HA 

9 O 

PA 

to 

1 

rnm 

O 

0 

1 

II 

o 

HA 

• 

Pm 1 

o’ 

| 

a 


m 

• IN 

- 
mm 

£ 

; s 

Hi 

0 sA i 

H P 

3 

I a 

1 

— , 
IN 

u 

s 

M 

X 

mm 

Nm 

N 

8 

M 

A 

PC 

5 

-N 

X 

M 

4D 

a 

m 

N* 

X 

43 

mm 

5 

IN 

i 

43 

T 

IN 

u 

s 

M 

A 

1 5 

H N 
43 A 
rnm »m 

"1 

IN N 

1 

f f : 

c 

N mm 

XX 

i : 

m *m 

C 

N 

f *° 

IN 

X 

M 

X 

"3f 

* 

: 

A 

5 

y 

: ■ 

i * 

A 

mm 

£ 

«m 

a A 

rnm 

'N. 

IN 

8 

: 

X 


1 

I d 4 
C Q u 

\ 2 ] 
* 4 i 

m U 
ft U. 

1 

£ 

• 

• 

4 

X 

am 

UJ 

• 

X 

1 

_£ 

1 

2 
A 

*m 

HI 

1 

_£ 

• 

£ 

i 

Jj 

• 

• 

2 

mrn 

£ 

2 

1 

4 

• • 

<§■ 3 

1 5 

3 4 

IJm-J • 

M 2 ! 

■ «l • 

l^ pj 3? i 

i I 

n ; 

* . 

• • 
m t>— 
ft k • 
. • u 

\§ 

• 

« 

£ 

Nm 

(al 

2 

_£ 

1 ‘ 

X , 

i 

4 

4 

2 

m- 

m- 

CmJ 


H 

s 

IB 

1 

, 

\ 

uO 

5 

1 : 

fo 

< 

1 

E 


h— 

..e 

3 

3 

1 

1 

stro^ 

« 

u 

i< 

r ^ 

r ; 

j r 

IN 

i 

M«n # o 

b 

4 

*3 — 

ft 

2 HA 
U Nr 

| 

« 

1 

F 

. 

IN 

y 




































































































































































































































160 


its peak. The displacements are given at eight stations along the 
blade. Also, the ratios of undamped and damped responses at all 
stations are given. This is one of the most important findings we 
have obtained so far, for it shows how much less the response of the 
damped system is as compared to that of the undamped system. For 
our models, at the excitation frequency of 400 [rad/sec], this 
ratio is about 1:?0. We cannot take this Information as the definite 
parameter of our system, it is only true when the system is driven 
at a frequency of 400 [rad/sec] . At different frequencies, this 
ratio will be different, and generally it will tend to Increase if 
we are driving the system closer to the natural frequency. For our 
models, this ratio is fairly high. Knowing that the low natural 
frequencies are closely spaced, when we are’ driving the system 
between the first and second natural frequencies we are fairly 
close to the resonance where damping becomes a very important parameter 
of the system. Also, for the thickest blade, we get a very uniform 
ratio between the undamped and damped responses at the stations 
along the blade. When thinner blades are used, this ratio is not 
so well (uniformly) distributed. Taking into account the 
geometrical and material characteristics of the blade, and the 
condition of the dynamic load (the amplitude of the excitation force 
is 1 lb.), we mieht be very close to the borderline of the elastic 
range of the blades. 

Also, another very important observation is that when the blades 
are driven (numerically) at a frequency less than the first natural 
frequency, the response of the blades at stations along the length is 
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very close to the static displacement of the load at that Instant 
of time. In other words, the response of the system at the frequencies 
below the firsit natural frequency is directly dependent on the 
stiffness of the system. At frequencies higher than the first but 
not close to the second natural frequency, the response depends on 
the mass of the system. This is very well documented in Table 8.6 
for all blade models. 

From the results above, we can make a very important statement 
about the resonant motion of mechanical vibratory systems with widely- 
spaced modes. Basically, each resonance is characterized by three 
zones of activity. At frequencies below the resonance, the activity 
of the structural system Is controlled by stiffness terms, which 
are proportional to the displacement and contain phases. At the 
resonant frequency, all restoring forces come from the damping 
terms, which are proportional to the velocity of the driving point. 

At frequencies above the resonance, the restoring force is 
dominated by the Inertial mass terms, which are proportional to the 
acceleration at the driving point. 


Table 8.6. Ratios between undaaped and daaped response when the forcing function reaches 
peak. The blade was forced by the F ■ slnu t at the tip, u ■ 400 (rad/sec ) . 
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Chapter 9 
CONCLUSIONS 

An analytical model and the experimental results have been 
obtained for prediction of the torsional stiffness of five blade 
models having NACA four-digit symmetrical airfoil cross sections. 

These five models, three with the airfoil cross sections and two 
with rectangular cross sections (aluminum and steel) , were modeled 
by modal analysis as either undamped or damped multi-degree-of- 
freedom systems. The complete identification of modal parameters 
for both cases was carried out. 

The dynamic response of undamped and damped blades was 
calculated. The excitation force F(t) * sin(w^t) was acting at 
the tip of the beam. The natural frequencies and lumping at 
particular modes were obtained by tasting with a Nicolet 660A Dual- 
Channel Vibration Analyzer. Based on the measured data, the 
analytical model was used to model damping in such a manner as to 
get a decoupled equation of motion. A small code was developed to 
identify all modal parameters and to predict the dynamic response 
of undamped and damped blade models which were modeled as multi-degree- 
of-freedom systems. 

From the analytical, experimental, and numerical phases of this 
study, the following conclusions can be made: 

1. Fairly accurate results were obtained for the torsional 
stiffness of blade models having NACA four -digit symmetrical airfoil 


cross sections. 
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2. The results for torsional stiffness have been scaled in 
non-dimensional form, such that for any model of an arbitrary 
length having a NACA four-digit cross section, the torsional constant 
can be calculated as product of the torsional stiffness, the scaling 
factor, and the reciprocal value of the length. 

3. Predicted values of natural frequencies were in good 
agreement with experimental results obtained by the FFT analyzer. 

4. The modal parameters were identified when either undamped 
or damped blade mode,ls were considered. They were essential for 
dynamic analysis of multi-degree-of -freedom systems. 

5. A good agreement was obtained between the continuous 
and lumped parameter models for the dynamic response of the blade 
when the tip of the blade was subjected to a harmonic force of 
F(t) ■ sin(w^t) . 

6. The dynamic response of the five blade models was predicted 
using the same excitation force. The damping obtained in the 
experimental phase was modeled analytically and Included in the 
calculation of dynamic response. 

7. A comparison of undamped and damped responses of the 
blades «as mad <T for the same load conditions. The ratio was 
determined, and was very much dependent on the forcing frequency 
and the amount of damping in the system. 

8. A much better insight into the damping of aluminum blades 
and its influence on the dynamic response of blades was obtained. 

This is of special interest when considering a system driven on 
high frequencies which are fairly close to resonance. 
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Several recommendations for future work can be made from this 
study, as follows: 

1. Since the torsional vibration mode has to be considered 
in propeller, rotor, and turbomachinery blade dynamics, accurate 
calculation of torsional stiffness is of special interest. This 
study could be used as a basis for the development of a small 
finite-element code using Galerkin's method. Experience has shown 
that this method is more accurate and universal than other approximate 
methods in boundary-value problems. References such as Richards (62), 
Vemuri and Karplus (88), Segerlind (73), and Rao (61) give a very 
elaborate treatment of boundary-value problems by employing the 
finite-element approach and Galerkin's approach. 

2. With the information on the amount of damping in the 
system provided in this study, this model can be related to another 
model which may be more convenient in application. The best way 

to relate the two models is to equate the dissipation of energy 
per cycle for both models. Knowing the damped dynamic response, the 
damping constants of dashpots located at the nodal points could be 
calculated . 

3. For the experimental phase, it is recommended that a 
modern FFT analyzer be adapted for modal analysis in (Structural 
dynamics, with an option in the methods for identification of 
modal parameters. Using these analyzers, one would be able to 
visualize the system and to observe the deformation pattern during 
different load conditions (transient and steady state). 

4. The concepts of this study for the analytical and 
experimental phases can be used for identification of modal 
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parameters in the torsional mode of vibration. With a slight 
modification, the computer program is capable of handling the 
torsional vibration on a fall scale, 

5, The most difficult task one faces is to develop gradually 
the analytical and computational procedure for dynamic response of 
real blades. In order to accomplish that, one must consider the 
complicated geometry of blades, the strong coupling between 
flexural and torsional nodes of vibration, centrifugal stiffening, 
the effect of the Coriolis force for thin and long blades, and 
unsteady excitation forces varying along the length and azimuth 


of blades 
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Appendix A 

MODAL PARAMETERS IDENTIFICATION 
A.l. Introduction 

Before we start considering the techniques for modal parameters 
identification, it is instructive to describe the distribution of 
the magnitudes of the transfer functions along the frequency range 
in real structural systems. 

When a structure is excited by a broad-band input force, many 
of its modes of vibration (degree, '•-of-freedom) are excited simul- 
taneously. Since the structure is assumed to behave in a linear 
manner, its transfer functions are really the sums of the resonance 
curves for each of its modes of vibration, as shown in Figure A.l. 

In other words, at any given frequency, the transfer function 
represents the sum of all the modes of motion which have been 
excited. However, near the natural frequency of a particular mode, 
its contribution to the overall motion is generally the greatest. 

The degree of mode overlap, i.e,, the contribution of the tails of 
adjacent modal resonant curves to the transfer-function magnitude 
at a mode's natural frequency, is governed by the amount of damping 
of the modes and their frequency separation. Figure A.l illustrates 
light and heavy modal overlap. Figure A.l(b). shows modes with 
light damping and sufficient separation so that there is a little 
modal overlap. Figure A. 1(c) shows modes with heavy damping 
such that there is plenty of modal overlap. 

In cases where modal overlap is light, the transfer function 
data can bq considered in the vicinity of each modal resonance as if 


Magnitude 
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it were a single-degree-of -freedom system. In other words , it is 
assumed that contribution of the tails of adjacent modes near each 
modal resonance is negligibly small. On the contrary, when modal 
overlap is heavy, a single-degree-of-freedom for modal parameter 
identification will not work; the parameters of all the modes musf 
be identified simultaneously. 

Usually, four residue-extraction techniques are used and named 
by the method: 

1. Quadrature response 

2. Peak picking 

3. Circle fitting 

4. Least squares. 

In the following sections, a detailed treatment and application is 
given for all four medal parameter Identification techniques. 


A. 2. Mode Identification by the Frequency Response 
The accuracy of the residue estimate obtained by this method 
is a function of how well each mode of the system is represented 
by a single-degree-of-freedom model. The magnitude of the complex 
frequency function for the r-th mode |h(.w ) | is determined from 

|h(jw r )|- / h(03 r ) -h*(io r ) CA.l) 


where h (jW r ). is the complex conjugate of the frequency response 
function 


h*(jw r ) 


1 ' <sf> - 

r r 


(A. 2) 
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o* Qimm, 

Thus , , 

R ^ 

I 2 — r 2~ (A. 3) 

+ < 2 Vr» 

where R is the modal residue and ? and ft are defined by the 
. t r r 

control theory approach. The plot of |h(w r )| is shown in Figure A. 2, 
The peak amplitude of the frequency-response function is obtained 
by differentiating eq, (A. 3) with respect to W and letting the result 
equal zero. The peak occurs when 

W r - ft r / 1-2 ? r i . (A. 4) 

Thus, the peak of the magnitude of the frequency response function 

occurs at a frequency less than both U)^ and ft^, For very lightly 

damped systems, the peak occurs at approximately the undamped 

natural frequency ft . 

r 

The maximum value of |h(to r ) | for a single-degree-of-freedom 
system in the r-th mode is referred to as the "Q r " of the system. 

When a system has a light damping, 

l h( ‘V I * «r “ ?r- (A - 5) 

An approximate value of damping for a single-degree -of -freedom 
system can be computed by using the frequency difference between 
the system's half -power points. The half-power points are defined 
as the amplitude 


2 


Frequency Ratio 


Figure A. 2. Magnitude of a frequency-response function for a 
singie-degree-of -freedom system (in the r-th mode) 
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or in terms of the magnitude |h(w)|, when it is equal to "Q^" . The 
half -power points occur when 


Q r 2 

<rj* 

ST 


2 V 2C ' 
r 


(1 - (£> V + ( 2 ?r ^.) 2 

r r 


(A. 7) 


Assuming light damping in the r-th mode, C r « 1. eq. (A. 7) becomes 

2 


(i) 


(#-) -1 + 2V 


(A. 8) 


Now, let the roots of eq. (A. 8) be equal to the frequencies u>^ and 
b>2« Then the difference of these two frequencies is referred to 
as the bandwidth of the resonance 


OJ 2 - ■ Aw. - bandwidth. 


For lightly-damped structures. 


Au> * 2 5 ft . (A. 9) 

r r r ' 

The approximate damping ratio of the system can be obtained from 
the magnitude of the frequency-response function by 

Aw 



(A. 10) 


r 


The identical procedure can be repeated for the other modes. We 
derived the formulas above under the condition of light modal 
overlap (the modes are well spaced along the frequency range) . 
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A. 3. Mode Identification by the Peak Picking Method 
The peak picking method is essentially a single-degree-of- 
freedom system technique. Again, as in the previous case, the 
accuracy of the residue estimate is a function of how well each 
mode of the system is presented by a single-degree-of -freedom 
system model. From system dynamics it is known that the properties 
of the transfer function for a single-degree-of-freedom system in 
the r-th mode could be determined from measurements of the frequency 
response function, expressed in terms of its poles and residues as: 

* 


K R 

j* 

h(t V r! 2j (jw r -p r ) " 2jTjur=pjr- 


(A. 11) 


For the single-degree-of-freedom system, the residue is seen to be 
real value, given by 


1/m 


R_ - - 

r 0) / 1 - C 
r r 


•mtx) 


(A. 12) 


dr 


In general, the residue may be considered complex and may be 
expressed as 


\ " l\|e r - |R r |(a r +jb r ), 


(A. 13) 


where j is the magnitude of the residue 


|R r l - / Re 2 (R r ) Im 2 (R r ) 


(A. 14) 


and a is the phase angle of the residue in the r-th mode 

Im(R ) 

“r ‘ arctan IS(rT • 


(A. 15) 
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The constants a and b are expressed as 


a 


r 


b 


r 


* COSd 

r 

■ sintt^. 


(A. 16) 


Substituting the expressions above for R r in eq. (A. 11), we have 


h(U) f ) 


l R J 

a (to. -to)+b C 0) 
r dr r r r 

|r 1 

r 

b (to, -to) -a ? to 
r dr r r r 
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a (to, 4to)+b ? to 
r dr ' r r r 

3 2 

b (to. -Ho) -a ? to 
r dr ' r r r 

2 

(to dr -Ko) 2 +(; r to) 2 

(to dr -HO) 2 +(C r to) 2 


(A. 17) 


In eq. (A.17), the real and Imaginary parts of the frequency-response 
function are explicitly revealed. Again, the first two terms 
represent the contribution of the positive pole, while the last two 
terms represent the contribution of the negative pole. 

For the purpose of graphic interpretation, it is useful to 
rewrite equation (A.17) in a dimensionless fora. To accomplish this, 
we divide eq. (A.17) by the static frequency-response function: 


h (0) = 
r 



/l-C 2 + 


? r“ + Mr 

r r r 


U) 



(A. 18) 


For a single-degree-of -freedom system, the static frequency-response 
function is simply the inverse of the spring constant. The following 
dimensionless frequencies are also defined as: 
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and 


e 




Equation (A. 17) may be rewritten as 


h(w r ) 

1 

q r (0-Y)+b r ? r 
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b r (8fy)+a r C r 

(&fY) 2 +C r 2 

" 2A 

(f*Y) 2 +C r 2 

, ,. TO 


— — — 


(A.19) 


(A. 20) 


Equation (A. 20) is shown graphically for various values of the damping 
ratio in Fig. A. 3. For all curves in this figure, the residue is 
taken as real (a " 1, b ■ 0). Thus, Fig. A. 3 represents the frequency- 
response function for a single-degree-of -freedom system. Figure A. 3(a) 
shows the real part of eq. (A« 20) , while Fig. A.3(.b) shows the imaginary 
part. It should be noted that the real part is symmetric in frequency, 
while the imaginary part is anti-symmetric. This property is called 
conjugate evenness, and is characteristic of structures with real 
residues. Since, for structural systems, we have interest only in 
a positive frequency range, eq. (A. 20) is plotted in Fig. A. 4 for 
positive frequencies only. In this figure, the residue is real and 
the damping ratio of 0.05 is selected. The parameters of interest 
in this figure are the frequencies at which the real and Imaginary 
parts of the frequency response function have extremes. 


Transfer Function (Imag. Part) Im [h(w )] Transfer Function (Real Part) Re[h(o) )] 
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- 2.0 - 1.5 - 1.0 - 0.5 0 0.5 1.0 1.5 2.0 


Frequency Ratio (to/u) ) 
(a) r 



Figure A. 3. Single-degree-of -freedom system transfer function (real 
modes); (a) real part, (b) imaginary part. 
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The frequency corresponding to the negative peak in Fig. A. 4(b) 
(imaginary part) is found by equating the frequency derivative of 
the imaginary part of eq, (A. 17) to zero. For a real residue, the 
imaginary part of eq. (A, 17) is 


I [h(o> ) ] - 
m r 


R r ' 

-S 0) 
r r 

R 

r 

w 

r r 

2 , 

(o) dr -o)) 2 +(r >r O)r ) 2 

" 2 

(0) dr +o)) 2 +(C r 0) r ) 2 


. (A. 21) 


Equating the derivative to zero gives 


_d_ 

da) 


I [h(w ) 1 
iii r 


0 



WV 0 



|R r UrVr (0)) • 


(A. 22) 


In eq. (A. 22), the term e^oi) represents the contribution of the 
negative pole. For frequencies near the natural frequency of the 
system, this term is negligible 


e (U*w ) = — ^ . (A. 23) 

r r 80) 3 
r 

Neglecting this term, equation (A. 22) indicates that the peak in 
the imaginary display of Fig. A. 4(b) occurs at 

0) * 0), . ’ (A. 24) 

dr 


That is, the peak in the imaginary part of the frequency response 
function is at the damped natural frequency of the system. If we 
consider the real part of the frequency response shown in Fig. 

A. 3(a), near the damped natural frequency this function has two 
extreme values, one positive and one negative. The frequency at 
which these extremes occur is found by equating the frequency derivative 
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of the real part of eq, (A. 17) to zero. For a real residue, the 
real part of eq. (A. 17) is 


R e th(u r )] 


R 

r 

(to -w) 
dr 

R 

, JP 


2 

(a) dr -w) 2 +(s r w r ) 2 

+ T 

(0) dr -H0) 2 +(? r (0 r ) 2 


(A. 25) 


Equating the derivative ‘of this expression to zero gives 


— 6 r (w) , (A. 26) 

In the expression above, • the term e^W) again represents the contribu- 
tion of the negative pole. For frequencies near the natural frequency 
of the system, this term is negligible: 

e (wsu) ) s -L-, (A. 27) 

r r 4uf 

Neglecting this term, eq. (A. 26) Indicates that the extremes in the 
real display of Fig. A. 4 (a) occur at 

w - “ dr + ? r U>. (A. 28) 

The difference between these two frequencies is 

Aw« £ to . 
r r 


XT R th(co ) ] - o - -~ 
dco e r 2 


(u dr -“) 


(u -w 2 )+(; u ) 2 

dr r r 


The damping is then estimated from 



(A. 29) 
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The damping properties at resonance In the r-th mode can be described 
by an amplification factor Q^. The amplification factor is 
related to the viscous damping c"or, C, by the relationship 


Q 


r 



(A. 30} 


where is the damping ratio of actual and critical damping of the 
system in the r-th mode. This method is the most precise in practical 
Implementation, particularly when dealing with high Q modes. It 
has the distinct advantage of being directly applicable to the real 
part of the transfer function. From Fig. A. 4 (a) , we have 


w ; + 1 

r (u /w ) 2 - 1 

ar br 


(A. 31) 


where 

0)^ - frequency above resonance in the r-th mode, where the real 
part of the frequency transfer function reaches the peak. 

■ frequency below resonance in the r-th mode, where the real 
part of frequency transfer function reaches the peak of the 
opposite sign. 

To summarize, the essence of the peak-picking identification method 
is, therefore, the determination of the three frequencies associated 
with the extreme in the real and imaginary representation of the 
frequency response function. This method is effective under the 
following conditions: 

1. The frequency response function has real, or very 


nearly real, residues. 
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2, The frequency response function can be approximated by the 
sum of the frequency response functions of several single-degree-of- 
freedom systems whose natural frequencies are distinct and "well 
separated." 

3. The Influence of the negative poles on the measured 
(positive frequency) frequency response function is negligible, for 
moderate and high damping may not be consistent with the accurate 
determination of the modal damping ratios. 

Whether or not adequate separation of poles exists in a given 
frequency response function is a matter of judgment and experience. 

In general, a minimum of five half-power bandwidths (Aw) should 
separate modes for adequate measurement purposes. Two extreme 
examples are shown in Figs. A. 5 and A. 6. 

Figure A, 5 shows two "well-separated" modes. Between resonances, 
both the real and imaginary parts of this frequency response function 
are relatively flat. Figure A. 6 shows two "poorly separated" modes. 
The effect of one mode on the other is quite apparent, especially in 
the real part of the frequency resp >tse function. Similar phenomena 
can be observed if the magnitude of the frequency response function 
is plotted. 



Transfer Function (Imag. Part) Im [h(u) ] Transfer Function (Real Part ) Relh(to) ] 
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(b) 


Figure A. 5. Tvo-degree~of -freedom system transfer functions with two 

distantly coupled modes; (a) real part, (b) imaginary part 



Transfer Function (Imag. Part) Im [h(to) ] Transfer Function (Real Part) Re (h(w) ] 
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(a) 



Figure A. 6. Two -degree-of -freedom system transfer functions with two 
closely coupled modes; (a) real part, (b) imaginary part. 
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Appendix B 

SYSTEM EXCITATION AND EXPERIMENTAL TECHNIQUE 
B.l. Introduction 

The transfer characteristics may be measured by isolating the 
structure from its normal environmental loading and subjecting it 
to an excitation of a simpler nature. In practice, the structure 
is subjected to a force introduced at a single point with a fixed 
orientation, and the resultant motions ore measured at one or more 
fixed points. This permits measurement of a series of transfer 
functions, thus characterizing the structure. 

Each transfer function provides two pieces of cause/effect 
information: 

1. The ratio of force/motion (as a function of frequency) 
relating the force input point and the motional response 
point* 

2. The degree of phase lag (as a function of frequency) 
between the force and motion introduced by the structural 
path. 

As we stated previously, the transfer functions identify the 
resonant frequencies, damping characteristics, mode shapes, stiffness, 
and inertia of the structure. 

Six transfer -function types are frequently measured. These are: 

1. Compliance = Displacement/Force 

2. Mobility = Velocity/Force 
Inertance = Acceleration/Force 


3. 
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4* Apparent Stiffness - Force/Displacement 

5. Impedance ■ Force/Velocity 

6. Apparent Mass ■ Force/Acceleration 

All six of the descriptive transfer functions are interrelated by 
simple algebraic operations. Having knowledge of one of the six 
frequency response forms is equivalent to having knowledge of all 
six. The natur.e of available instrumentation, local conventions 
and history, and problem specifics normally dictate the selection 
of one of the six functions for measurement. 

B.2. System Excitation 

The transfer functions can be measured using three broadly 
different test procedures. These techniques are called swept sine, 
random, and Impulsive excitation. There is no "right" technique 
for all problems. Modern FFT analyzers are fully capable of 
performing analyses for all three types of testing. The following 
is a summary of all three techniques. 

Swept-sine testing is the oldest of the three techniques. In 
this procedure, the structure is subjected to sinusoidal excitation, 
and the frequency of the sinusoid is slowly swept through a range 
of frequencies. If the sweep race is sufficiently slow (as 
determined both by the structure and the measurement equipment) , 
the transfer functions computed from this technique are identical 
to the results obtained by "resting" at a sequence of frequencies. 

i 

In random excitation, testing for the source of excitation is 
provided by a broad-band noise generator which produces an output 
signal with a "white" spectrum (a spectrum flat to at least the 



192 


bandwidth of the analysis)* Probability distribution of the broad- 
band excitation is Gaussian. This source of excitation is fed 

through a mass-compensation network, if required, for the examination 
* 

of the light structures with high "Q" modes, then through a power 
amplifier to the drive coil of the shaker. Both the excitation 
force and the response variable are measured simultaneously using 
a dual-channel FFT analyzer. 

The methodology for impact testing involves striking the 
structure with an impactor, usually a hand-held hammer, which has 
a load cell attached to measure the force input. The input and 
output, usually accelerations, are simultaneously measured. The 
two signals are Fourier-transformed, and the frequency response 
ratio is formed. It is assumed that the impulse is of sufficiently 
short duration to have the necessary frequency content to excite 
all frequencies in the band of interest simultaneously. Since this 
technique is very much of interest in this area of study, it is 
going to be treated in detail in the following sections. 

To summarize, sinusoidal testing is the slowest but most 
precise procedure for performing frequency response measurements. 
Impulsive force-input techniques are the fastest but least precise. 

Sinusoidal techniques have retained their popularity throughout 
history. Transient techniques have come into recent vogue as a 
means of performing rapid analysis on complex problems. The quick 
"look-see" capability offered by impulsive techniques is often 
sufficient to solve a problem. Frequently, initial studies are 
conducted using impulsive techniques, and higher -resolution information 
is gathered (where required) using either the random or swept-sine 
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techniques, 

FFT processing analyzers have gained enormous popularity because 
they are capable of performing transfer-function measurements using 
any of the three types of excitation. 

B.3. Impulsive Testing 

This procedure utilizes short-duration Inputs (with corresponding 
broad-band spectra) to excite all frequencies in the structure 
simultaneously. Impulsive testing is normally conducted using 
some form of an instrumented hammer such as that shown in Fig, B.l. 

A force transducer is mounted either on the head of the hammer or 
on that portion of the structure to be impacted. 

An accelerometer is used to measure the response of the 
structure. By far the most popular technique for impulsive 
excitation is to mount an accelerometer at a fixed location and to 
excite the structure at a multiplicity of locations using a hammer 
with a force transducer on its face to impact the system at many 
points. The data is analyzed with a dual-channel FFT analyzer. 

Each time the structure is impacted with the hammer, both the 
impulse and the response acceleration are captured by the FFT 
analyzer, using the impulse force as a transient-capture trigger 
condition. These t.vo transient histories are Fourier-transformed to 
yield the input a nj output spectra. The resultant ratio of the 
output and input spectra is the desired transfer function. 

Usually, the results of several transient excitations are 
averaged. In this situation, the most desirable procedure is to 
compute the cross-spectrum between the input and response, and the 
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Figure B.l. The setup for Impulsive testing. 
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power spectrum of the input. The resultant averaged transfer 
characteristic is derived by dividing the cross-spectrum by the input 
power spectrum. 

There are several advantages to impulsive testing. It is, 
by far, the most rapid technique that can be applied to the study of 
a complex structure. Its fundamental advantage is that no shaker 
system needs to be mounted to the test object. The testing time 
required for each spacial position on the structure is determined 
by the number of times the tester wishes to impact the structure 
with a hammer. Because the response transducer is mounted at a 
fixed location, no time is required to move the response transducer 
from point to point and, of course, a hammer may be easily moved 
from one point on the structure to another. 

Some disadvantages of Impulsive testing, however, also exist, 
as follows: 

1. Impulsive techniques impose stringent requirements on the 
analysis instrumentation. Specifically, the FFT analyzer must 
have considerable more dynamic range than is required for sinusoidal 
or random testing. 

2. Precision of the results is strongly affected by the care 
with which data is initially acquired. If the signal-conditioning 
amplifier or analog-to-digital converters of the FFT analyzer are 
allowed to momentarily overload during the acquisition of a 
transient, the results are grossly affected, 

3. It is important to match the duration of the impulsive 
force input to the band of frequencies to be analyzed. This requires 
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selecting the mass and geometry of the hammer as well as the nature 
of the material and the hammer /structure interface to "shane" 
the Input pulse. 

4. The mass of the striking "tip" between the force transducer 
(mounted on the face of the hammer) and the structure can introduce 
errors when exciting light structures. Unfortunately, there is no 
currently-available technique for performing on-line mass cancellation. 
The results of impulsive testing have lower signal-to-noise ratios 
than those from any other form of excitation. 

B.4. Impact or Selection 

The Impactor is usually a hand-held hammer instrumented with 
a force transducer. The force transducer has an impact cap attached 
to it, and the hammer may have an extension attached to it. The 
hammer could be very small, weighing only a few ounces, or as large 
as a sledge hammer. 

The frequency content of a force pulse is determined by its 
shape. Often it is desirable to shape the pulse to ensure that 
the impact excites the structure in a particular desired frequency 
range. The three parameters which characterize the pulse are: 

1. Pulse height . This is largely controlled by the impact 
velocity. A higher impact velocity would require a larger force to 
decelerate a given mass, resulting in an increase in pulse height. 

2, Rise time . This can be controlled by the choise of impact 
caps. If the impact cap is thought of as a spring, it seems 
reasonable that an impact with a soft spring (soft impact cap, such 
as rubber or plastic) would have a longer rise time than an impact 
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with a hard spring (hard impact cap) . 

3. Pulse width . To a certain extent* the pulse width can be 
controlled by the use of extenders to change the hammer mass. For 
example, a larger mass (hammer with extender) moving at the same 
velocity will have a larger amount of kinetic energy and require a 
longer time to be stopped by a given spring (the Impact cap) , 

This produces a wider pulse. 

It is necessary to point out that methods of shaping the pulse 
are all Interactive, For example, installing a softer impact cap will 
increase the pulse width as well a;J Increase the rise time. It is 
best to try various configurations of impact caps, extenders, and 
Impact velocities and Fourier-transf orm the time signals to ensure 
that the desired frequencies are being excited. 

B.5. Hammer Calibration 

In order to calibrate the frequency response plots generated 
by impact testing, it is necessary to dynamically calibrate the 
load cell which is part of the hammer structure. In Fig, B.2, a 
typical calibration setup for an impact hammer and accelerometer is 
shown . 

On impact, the structure experiences the force which is 
necessary to stop the entire hammer assembly. The load cell’s 
crystal element experiences only the force necessary to stop the 
mass of the hammer head and extender. It does not sense the force 
component necessary to stop the mass of the impact cap and the 
load cell's seismic mass, Because of this, the sensitivity of the 
load cell appears to change as the mass of the impact cap and the 
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Figure B.2(a). Typical calibration setup for an impact hammer and 
accelerometer set, (b) typical impact hammer with 
a load cell. 
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mass of the hammer are changed. Increasing the ratio of the hammer 
head plus extender to the total mass of the hammer structure 
decreases the apparent sensitivity change, but as the total mass 
increases, it can introduce problems of multiple impacts and 
penetration. The mathematical model of hammer calibration is given 
in the following text. 

The hammer calibration involves testing of the functional 
transfer behavior of a test object, with a known mass acting as a 
rigid body employed as the test object. For scaling a test object 
measurement, only the ratio of the accelerometer sensitivity to 
the hammer sensitivity (S /S.) is needed. However, the actual 
sensitivity of Che hammer "S^" can be determined from equation (b) 
in Fig. B.3 when a precision quartz accelerometer with a known 
sensitivity "S " is employed as a reference. The "sensitivity" of 

a 

a load cell or hammer structure is just another way of expressing 
its functional transfer behavior. Calibration (scaling) factors 
can be calculated from the peaks of the time-varying signals from 
the ratio of their frequency components. 

The sensitivity of a hammer differs from the static sensitivity 
of its force-sensor component because of the normal behavior of the 
hammer structure. The force experienced by the crystal sensing 
elements sandwiched within the hammer structure is not the same 
as the force applied on the test object. For this reason, the hammer 
assembly ought to be calibrated as it is used. To test small 
hammers, a several-pound cylindrical mass with the accelerometer 
attached and suspended with a long leads should be used. For a given 
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hammer structure, the sensitivity of the moving hammer differs from 
its static sensitivity by a constant percentage proportional to 
the ratio of the mass of the tip and cap in front of the crystals to 
the total mass of the hammer assembly, 

Three methods of scaling a frequency response measurement are: 

1. Inserting the average sensitivities supplied by the 
manufacturer. 

2. Measuring and using the "ratio'* of sensitivities as out- 
lined above, 

3. Calibrating the hammer sensitivity, using an accelerometer 
as a reference. 

With methods 2 and 3, an average scaling factor can be estimated 
or a scaling factor at each discrete frequency in the spectrum of 
interest can be measured and used, thus compensating for any 
anomalies in the hammer or motion-sensor behavior. 

Behavior objectives for hammer and sensor structures are straight 
lines relating input and output and their ratio (sensitivity) with 
frequency. In other words, ideal sensors treat amplitudes proportion- 
ally and frequencies the same, and do not delay the signal. The 
results for one typical Impact hammer are given in Fig. B.4. 
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Forc« 

(c) 


Figure B.4. Typical results from calibration of an impact hammer and 
accelerometer; (a) time record, (b) frequency response, 
(.c) amplitude linearity. 
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Appendix C 


FLOW CHART FOR THE COMPUTER PROGRAM FOR PREDICTION OF 
DYNAMIC FLEXURAL RESPONSES OF PROPELLER BLADES 


Start 


Main program 
Call subroutine INPUT 


Subroutine INPUT 

1. Read and print material and geometrical properties of model, 

2. Read and print the amplitudes, frequency, and time Increment 
of the forcing function, 

3. Read and print the Initial displacements and initial velocities 
at nodal points. 

4. Calculate the cross-sectional area and cross-sectional moment 
of inertia for an N.A.C.A. four-digit symmetrical airfoil 
sectional area. 


Main Program 

Calculate and print lumped masses at all nodal points. 

Calculate and print flexibility Influence coefficients. 

Call function SIMUL. 

Calculate and print the stiffness influence coefficients matrix 
by inversion of the flexibility influence coefficients matrix. 


Call function SIMUL 


Function SIMUL 

When INDIC is negative, SIMUL computes the inverse of the N by 
N matrix in plate . 

When INDIC is zero, SIMUL computes the N solutions X1(1)..X(N) 
corresponding to the set of linear equations. 

When INDIC is positive, the set of equations is solved, but the 
inverse is not computed in place. 


m 
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Subroutine JACOBI 

This subroutine solves the generalized eigenproblem using the 
generalized Jacobi iteration. 

Input variables 

A(N,N) ■ Stiffness matrix (assumed positive-definite), 

B(N, N) m Mass matrix (assumed positive-definite) 

X(.N,N) ■ Matrix storing eigenvalues on solution exit 

EIGV(N) - Vector storing eigenvalues on solution exit 

D(.N) - Working vector 

N - Order of matrices A and B 

RTOL ■ Convergence tolerance (usually set to 10**-12) 

NSMAX m Maximum number of sweeps allowed (set to 15) 

IFPR ■ Flag for printing during iteration 

EQ»0 No printing 

EQ*1 Intermediate results are printed 
IOUT ■ Output device number 

Output variables 

All eigenvectors are stored columnwise 
A(N,N) - Diagonalized stiffness matrix 

B(N,N) ■ Diagonalized mass matrix 

X(N,N) b Eigenvectors 

XS(N,N) - Scaled eigenvectors XS(I, J)TR*M1(I,J)*XS(I,J)-(I) 
XL(N,N) ■ Orthonormalized eigenvectors 
EIGV(N) * Eigenvectors 


The eigenvalues and corresponding eigenvectors are not calculated 
in an ascending or descending order. 


Main program 


Call subroutine MODAL 
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Subroutine MODAL 

Tills subroutine defines and prints out the following parameters; 

1. The coefficients of the modal equation. 

2. The orthogonality relations of eigenvectors. 

3. Modal mass and modal stiffness. 

4. Effective mass in a mode. 

5. The percentage of total mass taken into consideration by modal 
analysis. 

6. The participation factor (modal force). 

7. Natural frequencies sorted in an ascending order, 

8. Critical damping, amplification factor, damping ratio, and 
damping coefficients in a particular mode. 


Main program 


Call subroutine SORT 


Subroutine SORT 

This subroutine defines and prints out the following parameters; 

1. The eigenvalues and corresponding eigenvectors sorted in an 
ascending order. 

2. Natural frequencies obtained from a continuous beam model. 

3. First five natural frequencies measured by FFT analyzer. 

Scaled eigenvectors sorted in an ascending order columnwise 
represent the modal shapes in a particular mode. 


Main program 


Call subroutine EXACT 
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Main program 

Call subroutine DRESP 

I 

Subroutine DRESP 

This subroutine calculates the dynamic response of a uniform 
continuous blade (cantilever beam) modeled by modal analysis, with 
the damping Included. The damping characteristics of the blade 
can be represented appropriately by using proportional damping in 
a mode superposition. An uncoupled system of equations Is used. 
For the complete response, the solution to all N equations, 

I * 1,2, ...,N, is obtained and the nodal point displacements are 
calculated by superposition of the response in each mode. Initial 
conditions are taken into account, and the transient response is 
added to the steady-state response. 

The test program calculates and prints out only the steady- 
state response of a blade subjected to a harmonic force of F ■ FO* 
SIN(OMEGAF*T) at the tip. Initial conditions are assumed to he 
zero (the system is driven from the rest). The displacements are 
calculated at N nodal points (along the length) at each time-step 
interval during two time periods of forcing frequency, 

♦ ; 


Main Program 




